{"id":"46c17e4e-4ff9-4610-8e4d-1d900aa197fd","arxiv_id":"2501.12856","paper_version":1,"verdict":"REJECT","confidence":"MODERATE","novelty_score":3.0,"correctness_risk":"high","formal_verification":"none","parameter_count":4,"one_line_summary":"The paper reformulates ODE parameter estimation as a derivative-matching least-squares problem and solves it with Gauss-Newton and gradient descent, showing better recovery than a bound-constrained NLS on three synthetic examples.","lead":"This paper estimates unknown parameters in systems of ODEs by matching finite-difference derivatives to the model equations, then applying Newton-Raphson and gradient descent updates. The authors claim these methods beat nonlinear least squares on test problems including a chaotic Lorenz system.","discovery_kind":"extension","skeptic_critique":{"model":"deepseek-v4-flash","headline":"The evidence for 'outperforming NLS' is not credible because the NLS baseline is constrained and appears stuck at bounds or local minima; a fair unconstrained multi-start NLS comparison could reverse the headline.","rationale":"I read the paper as claiming a practical advantage over NLS based on direct derivative matching. For that claim to hold, the comparisons in §3 must be fair and the derivative estimates must not be systematically corrupted by noise. The weakest link is the comparison protocol: NLS is constrained, no multi-start or solver details are given, and the reported NLS results are at or near bounds, so the 'outperformance' in Tables 3/7/10 may be an artifact of a poorly optimized baseline. A fair NLS fit to the same data would settle this directly. The reader's emphasis on finite-difference noise is related and valid; I mark partial agreement. I also noted a formal issue in Eq. (16): the update as printed multiplies an m×m matrix by an M×1 vector as (J^T J)^{-1} J E, which is dimensionally undefined unless m=M; the standard Gauss-Newton step would require J^T E. This supports rejection but is secondary to the unsupported comparison. The stochastic variants are never tested, and no code is provided, so the current paper should remain rejected; a fair-NLS check would determine whether the core numerical claim can be salvaged.","tokens_in":18450,"tokens_out":9667,"duration_ms":110349,"concrete_test":"Re-run the three §3 comparisons using a standard NLS solver (e.g., scipy.optimize.least_squares on the state-space sum of squared residuals) with unconstrained parameters and multiple initial guesses, including the true parameters and the paper's NR/GD estimates, on the same generated data and noise seeds. If NLS then recovers parameters closer to the true values or yields lower state RMSE than the reported NR/GD estimates, the claimed superiority over NLS is an artifact of the weak constrained baseline; if it still loses, the central claim survives this check.","verdict_should_be":"UNCHANGED","load_bearing_attack":"The central claim is that the proposed derivative-matching NR/GD estimators 'provide improved accuracy' and 'outperform NLS' (abstract). The only evidence for this is the three comparisons in §3.1.4, §3.2.4, and §3.3.4, and that evidence is not apples-to-apples: NLS is solved as a constrained optimization with parameter bounds (Eqs. (52), (57), (62)), while NR/GD are unconstrained. No multi-start, solver tolerance, or convergence diagnostics are reported for NLS, and the reported NLS solutions sit at or near the active bounds (a3=2 at the lower bound in Eq. (58); a3≈0.052 near 0 in Eq. (63)), indicating likely local minima rather than the actual NLS optimum. The error-metric tables (Tables 3, 7, 10) therefore do not establish that the proposed methods outperform NLS. Additionally, the proposed objective E in Eq. (14) is built from finite-difference derivatives (Eqs. (6)/(12)) of noise-added data; in the chaotic Lorenz case (§3.2) those derivative errors can be large, and the paper reports only single realizations, so the claimed accuracy is not demonstrated.","agreement_with_reader":"partial"},"referee_report":{"model":"deepseek-v4-flash","summary":"The paper proposes a direct derivative-matching approach for estimating parameters in systems of ODEs. The ODE residual is formed by comparing a finite-difference approximation of the state derivative (Eqs. (6) and (12)) with the model right-hand side f(t,x,a), and the parameters are updated by Newton-Raphson (Eq. (16)) or gradient descent (Eq. (30)) iteratively, with stochastic variants defined in Algorithms 2 and 4. The method is demonstrated on population dynamics, Lorenz, and activator-inhibitor examples, and the abstract claims that the proposed methods outperform nonlinear least squares (NLS). Convergence theorems are stated in Section 2.1 (Theorem 1) and Section 2.2 (Theorem 2).","tokens_in":18701,"tokens_out":3126,"duration_ms":31540,"significance":"If the central claim were established, the derivative-matching formulation would offer a computationally cheap alternative to NLS for ODE parameter estimation. However, the manuscript has load-bearing deficiencies: the convergence proof contains an impossible inequality, the numerical comparison to NLS is not apples-to-apples because the NLS baseline is constrained while the proposed methods are unconstrained, the stochastic variants are never numerically evaluated, and the treatment of noisy finite-difference derivatives is not analyzed. The paper does not currently provide a credible demonstration of its headline claim of improved accuracy over NLS.","major_comments":[{"comment":"Equation (39) states \"0 ≥ δ⊤δ = G(t, x, ai)⊤ Π² G(t, x, ai) + O(∥μi∥³)\". Since δ⊤δ is a sum of squares, it is always nonnegative, so the inequality cannot be correct as written; the only consistent statement would be 0 ≤ δ⊤δ. This sign error invalidates the derivation of the linear convergence bound in Eq. (41), which relies on the incorrect direction of the inequality. This is a load-bearing error in the GD convergence analysis.","section":"Section 2.2, Eq. (39)"},{"comment":"The claim that the proposed methods outperform NLS is not supported by the reported experiments because the comparison is not fair: the NLS baseline is solved as a constrained optimization with bounds (Eqs. (52), (57), (62)), while the proposed methods are unconstrained. The reported NLS solutions sit at or near the active bounds (a3=2 in Eq. (58) and a3≈0.052 in Eq. (63)), suggesting that the NLS solver is finding constrained local optima rather than the true unconstrained optimum. No multi-start, solver tolerance, or convergence diagnostics are reported for NLS. Therefore Tables 3, 7, and 10 do not establish the paper's central claim of outperforming NLS.","section":"Section 3.1.4, 3.2.4, and 3.3.4"},{"comment":"The stochastic versions SNR and SGD are defined in Algorithms 2 and 4, but the numerical examples only apply the deterministic NR and GD methods (Algorithms 1 and 3). No results are reported for the stochastic variants, despite the title, abstract, and conclusion emphasizing their contribution for large-scale systems. The statement in the introduction to Section 3 that the application of Algorithms 2 and 4 is presented in Python is not followed by any stochastic numerical experiment. Consequently, the claims about scalability and reduced computational cost for big data are unsupported.","section":"Section 3 (numerical implementation)"},{"comment":"The objective E in Eq. (14) is built from finite-difference derivatives of noisy data. The paper does not analyze how noise propagates through these derivative approximations or how that affects the estimated parameters. For the chaotic Lorenz system in Section 3.2, small state errors amplify exponentially, and the reported results are based on single realizations, so the accuracy claims are not demonstrable. Since the entire method rests on the reliability of these derivative approximations, this missing analysis is a load-bearing gap.","section":"Section 2, Eqs. (6) and (12)"}],"minor_comments":[{"comment":"The title contains a typo: \"ODE S\" should read \"ODEs\".","section":"Title and throughout"},{"comment":"The second equation writes \"dx2/dx\" instead of \"dx2/dt\".","section":"Section 3.1, Eq. (45)"},{"comment":"The step length γ is defined as the inverse of the Gram matrix of the Jacobian, which is appropriate only when the Jacobian has full column rank; the paper does not state this assumption.","section":"Section 2.1, Eq. (16)"},{"comment":"The citation \"Barzilai and Bowein, 1988\" should be \"Barzilai and Borwein, 1988\".","section":"Section 2.2, Eq. (31)"},{"comment":"The sentence \"The optimal estimated parameters ¯a int Table 4\" contains a typo: \"int\" should be \"in\".","section":"Section 3.2.3, text near Table 4"},{"comment":"The layout of the NR row in Table 4 is confusing because the convergent parameter vector is placed in the last cell, making it difficult to associate it with the corresponding initial guess.","section":"Section 3.2.3, Table 4"}],"recommendation":"reject","confidential_remarks":"The paper has serious mathematical and empirical problems: an impossible inequality in the convergence proof, an inequitable comparison to NLS, and unsupported claims about stochastic variants. Even though some issues could be addressed in a revision, the central claim of outperforming NLS is not supported by the present evidence, and the convergence analysis is not merely incomplete but contains a directly incorrect statement. The manuscript would require a substantial rewrite of both the theory and the numerical study to be publishable."},"author_rebuttal":null,"desk_editor":{"model":"deepseek-v4-flash","letter":"Short version: this is derivative matching, a setup that goes back to Varah (1982) and Ramsay et al. (2007). The paper's concrete additions are a finite-difference formula for non-uniform grids and the use of Newton-Raphson and gradient descent with a Barzilai-Borwein-style step size on that objective. That is a reasonable extension, but the evidence for the headline claim—'outperform NLS'—does not hold up.\n\nWhat it does well: the three numerical examples are worked through in detail, including a chaotic Lorenz system. Equation (12) for non-uniform differences is useful, and the GD step-size rule is sensible. The stochastic variants are clearly defined, even if they are never tested.\n\nThe soft spots are the usual ones, but they matter here. The NLS baseline is constrained with bounds while the proposed methods are unconstrained. In the Lorenz case the NLS solution sits at the lower bound (a3=2), and in the activator-inhibitor case a3≈0.052, so the NLS optimizer is likely stuck at a boundary or a local minimum. No multi-start, no convergence diagnostics, and all results are single realizations without error bars. The comparison tables therefore do not establish superiority.\n\nThe convergence analysis is sloppy. Equation (39) claims 0≥δ^Tδ, which is impossible for a nonzero vector; that is a sign error in a load-bearing inequality. The proof of Theorem 2 is not rigorous enough to support the stated linear convergence.\n\nThe literature review is broad and the writing is understandable, so the authors are not sloppy in intent. But the abstract overclaims. The stochastic versions are defined but never tested, so the 'big data' scalability statement is unsupported.\n\nWho is this for? Someone working on gradient matching could read it for the non-uniform grid formula and the chaotic-system test, but the method is not new and the comparison is not credible as is.\n\nMy recommendation: this deserves a serious referee rather than a desk reject, because there is a concrete method and testable claims. But the revision bar is high: fair unconstrained multi-start NLS comparisons, error bars, tests of the stochastic methods, and corrected proofs. Without those, the 'outperforming NLS' claim should not stand.","headline":"Derivative matching with a non-uniform grid twist, but the 'outperforms NLS' claim rests on a constrained baseline that is clearly stuck at bounds.","tokens_in":19258,"tokens_out":3355,"would_cite":false,"duration_ms":32550,"reading_group":"maybe","serious_thinker":"yes","would_accept_peer_review":true},"rs_alignment":null,"lean_confirmation":null,"pith_extraction":{"msc":["65L09","65K10","34A55"],"pacs":[],"model":"deepseek-v4-flash","headline":"This paper claims that ODE parameters can be estimated directly by fitting finite-difference derivatives to the system's right-hand side with Newton-Raphson or gradient descent updates, and that these direct methods beat nonlinear least…","keywords":["parameter estimation","systems of ODEs","Newton-Raphson","gradient descent","stochastic gradient descent","discrete derivative","Taylor expansion","nonlinear least squares"],"falsifier":"Generate Lorenz-system data with known parameters and additive noise, compute $x'(t)$ using the paper's equation (12), and check whether the NR/GD estimates converge to the true parameters as the number of data points $N$ grows. If the estimates remain biased even for large $N$, or if the residual objective's global minimum is not at the true parameter vector, the central claim fails.","tokens_in":18226,"feed_emoji":"📈","tokens_out":4736,"duration_ms":46301,"temperature":0.7,"pith_summary":"The paper tries to establish a simpler route to estimating unknown parameters in systems of ordinary differential equations: instead of repeatedly solving the ODE and matching predicted states to measurements, it replaces the derivative $x'(t)$ by a finite difference computed from the data and solves the residual equation $f(t,x,a)-x'(t)=0$ with Newton-Raphson or gradient descent. This direct formulation avoids numerical ODE integration during fitting, and the paper argues it yields parameters that are closer to the true values than those from the standard nonlinear least squares approach. The authors test the idea on a predator-prey system, the chaotic Lorenz system, and an activator-inhibitor model, reporting that Newton-Raphson converges in very few iterations while gradient descent handles chaotic dynamics robustly, and that stochastic subsampled versions scale to large data. If the claim holds, it offers a fast, simple alternative for fitting ODE models to noisy observations.","feed_headline":"Derivative fitting beats least squares for ODE parameters","feed_subtitle":"Newton-Raphson and gradient descent on a discrete-derivative residual recover parameters accurately, even for chaotic systems.","key_machinery":"The object that carries the argument is the residual vector $E(t,x,\\bar a)=f(t,x,\\bar a)-x'(t)$, with the derivative $x'(t)$ obtained from a finite-difference formula (first-order from two consecutive time points, or second-order from three points). Linearizing $f$ by Taylor expansion converts the residual into a linear system in $\\Delta a$, enabling the Newton-Raphson update $a_{i+1}=a_i-\\gamma\\,\\nabla_a F\\,E$ and the gradient descent update $a_{i+1}=a_i-\\eta\\,\\nabla_a G\\,G$. The machinery includes the Jacobian $\\nabla_a F$ (static for systems linear in parameters, updated each iteration otherwise), the explicit step-size rule for GD, and the random subsampling that defines the stochastic variants SNR and SGD.","core_discovery":"On the paper's own terms, the central discovery is that parameter estimation for systems of ODEs can be posed as a direct residual-minimization problem in parameter space using discrete derivatives and Taylor expansion, solvable by Newton-Raphson (NR) and gradient descent (GD). The residual is $E(t,x,\\bar a)=f(t,x,\\bar a)-x'(t)$, where $x'(t)$ is approximated by the first- or second-order finite differences in equations (6) and (12); expanding $f$ in a Taylor series in $a$ gives the linear system $E=-\\nabla_a F\\,\\Delta a$, which NR solves with the pseudo-inverse step $\\gamma=(\\nabla_a F^\\top\\nabla_a F)^{-1}$. GD instead uses the explicit step size $\\eta=\\delta^\\top G/\\delta^\\top\\delta$ with $\\delta=\\nabla_a G\\,\\nabla_a G^\\top G$. Stochastic versions randomly select $r$ of the $n$ ODE equations each iteration to cut cost. Across three numerical examples, the estimated parameters from both methods are closer to the generating parameters than those from constrained NLS, with error metrics closer to zero and $R^2$ closer to one.","pith_inferences":["Beyond the paper: the direct derivative-matching objective is a variant of gradient matching, and its bias under measurement noise is well known; a natural extension would be to derive the asymptotic bias of the NR/GD estimates and to add smoothing or regularization to the derivative estimate.","Beyond the paper: the authors compare against NLS with the same initial guesses, but do not report the basin of attraction for GD on the Lorenz system; measuring how often random restarts converge to the true parameters would sharpen the robustness claim.","Beyond the paper: the stochastic variants randomly subsample equations, not data points; an alternative extension is to subsample time points or use mini-batches of the residual, which could further reduce cost on very long time series.","Beyond the paper: the finite-difference derivative amplifies noise in chaotic systems, so a testable extension is to repeat the Lorenz experiment at several noise levels and check whether the parameter estimates remain closer to truth than NLS when the derivative is noisier."],"forward_implications":["Newton-Raphson on the direct derivative residual can converge in one or two iterations for systems linear in parameters, giving very fast estimates from clean or mildly noisy data.","Gradient descent with the explicit step size provides a robust fallback for chaotic systems, where Newton-Raphson may jump to a poor local solution.","Stochastic subsampling of the ODE equations makes the method scalable to large data sets by reducing the per-iteration Jacobian and residual size.","If the accuracy comparison to NLS holds, the direct method is a practical choice for real data fitting when the ODE lacks an explicit solution and numerical integration is expensive.","The finite-difference formulas accommodate irregular sampling intervals, so the method applies directly to experimental data with varying time steps."],"supporting_citations":[{"why":"Introduces the spline least-squares derivative-fitting approach that the direct method builds on.","marker":"[Varah, 1982]"},{"why":"Develops the generalized smoothing and profiled estimation alternative against which the direct method's simplicity is positioned.","marker":"[Ramsay et al., 2007]"},{"why":"Shows how multiple shooting fits ODEs to chaotic data, providing the comparison the paper's robust GD result addresses.","marker":"[Baake et al., 1992]"},{"why":"Supplies the explicit step-size update and descent guarantee used in the paper's gradient descent algorithm.","marker":"[Hao, 2021]"},{"why":"The two-point step-size gradient methods that motivate the choice of the GD step size.","marker":"[Barzilai and Bowein, 1988]"},{"why":"The multiple shooting method for ODE parameter estimation that the direct method is compared against in spirit.","marker":"[Peifer and Timmer, 2007]"}],"fun_headline_variants":["Discrete derivatives let Newton-Raphson and GD fit ODEs better","Stochastic Newton-Raphson and gradient descent sharpen ODE estimation","Newton-Raphson and gradient descent beat least squares for ODEs","Taylor-expanded residual lets NR and GD fit ODE parameters","Stochastic NR and GD improve ODE parameter estimation accuracy"],"cache_read_input_tokens":3200,"weakest_assumption_plain":"The finite-difference derivative computed from noisy data is close enough to the true derivative that the residual $f(t,x,a)-x'(t)$ has its minimum at the true parameters, an assumption that is especially fragile for chaotic systems where state errors amplify.","fun_headline_variants_meta":{"raw":{"variants":["Discrete derivatives let Newton-Raphson and GD fit ODEs better","Stochastic Newton-Raphson and gradient descent sharpen ODE estimation","Newton-Raphson and gradient descent beat least squares for ODEs","Taylor-expanded residual lets NR and GD fit ODE parameters","Stochastic NR and GD improve ODE parameter estimation accuracy"]},"model":"deepseek-v4-flash","effort":"low","cost_usd":0.000932,"raw_usage":{"total_tokens":4050,"prompt_tokens":1065,"completion_tokens":2985,"prompt_tokens_details":{"cached_tokens":384},"prompt_cache_hit_tokens":384,"prompt_cache_miss_tokens":681,"completion_tokens_details":{"reasoning_tokens":2894}},"tokens_in":681,"tokens_out":2985,"duration_ms":22748,"temperature":1.0,"reasoning_tokens":2894,"cache_read_input_tokens":384,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-08-10T16:44:05.869274+00:00","model_set":{"reader":"deepseek-v4-flash"},"falsifier":"Generate Lorenz-system data with known parameters and additive noise, compute $x'(t)$ using the paper's equation (12), and check whether the NR/GD estimates converge to the true parameters as the number of data points $N$ grows. If the estimates remain biased even for large $N$, or if the residual objective's global minimum is not at the true parameter vector, the central claim fails.","supporting_citations":[{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"The two-point step-size gradient methods that motivate the choice of the GD step size."}],"review_version":1}