{"id":"5e2321a1-a8af-4afe-afb6-2a2ea238992d","arxiv_id":"2411.15637","paper_version":3,"verdict":"CONDITIONAL","confidence":"MODERATE","novelty_score":5.0,"correctness_risk":"medium","formal_verification":"none","parameter_count":6,"one_line_summary":"A differentiable particle filter with L1 proximal updates estimates sparse polynomial transition functions and interaction graphs for nonlinear state-space models.","lead":"GraphGrad learns a sparse polynomial formula for the hidden dynamics of a noisy system from observations alone, using a differentiable particle filter and a penalty that forces most coefficients to zero. It recovers which components interact and how, on chaotic test systems like Lorenz 63, without ever seeing the true state directly.","discovery_kind":"new_method","skeptic_critique":{"model":"deepseek-v4-flash","headline":"The central recovery claim rests on treating stop-gradient DPF gradients as valid descent directions, but no bias/variance analysis or K-sensitivity evidence is provided; Eq. (15) may be minimizing a different objective.","rationale":"The reader's weakest assumption is that the optimization actually finds the penalized MLE: the stop-gradient differentiable particle filter is treated as providing a usable descent direction, and the proximal-gradient scheme is assumed to converge on a nonconvex, biased-gradient problem. This is exactly the most load-bearing point for the paper's central claim. If the gradient estimate is biased or the optimizer stalls in a poor local minimum, then the recovered C and adjacency graph A do not reflect the data-generating system even with unlimited observations. The paper provides strong empirical evidence on Lorenz 63, Lorenz 96, and Kuramoto, and the batched B-GraphGrad is a reasonable engineering solution to likelihood degeneracy. The comparison against the dense pMLE baseline is informative, and the Kuramoto model-mismatch experiment is a good stress test. However, no code is supplied and the inner-loop step count S is not reported, so the empirical claims cannot be independently reproduced. The missing bias/variance analysis for the DPF gradient is not an internal inconsistency, but it is a correctness risk: the optimization part of the pipeline is the least secure link between the method and the claimed inference of structure and values. A concrete K-sensitivity and gradient-angle check would settle whether the reported recovery is robust or an artifact of a particular particle budget. This concern does not overturn the paper; it makes the current CONDITIONAL verdict appropriate, so no change to the reader's verdict is recommended.","tokens_in":24940,"tokens_out":6761,"duration_ms":72300,"concrete_test":"Re-run the Lorenz 96 experiment of Sec. V-C with K = 500 and K = 2000 particles (keeping T, S, B, λ fixed) and compare the F1 and RMSE distributions over the same 150 seeds. In parallel, at a fixed C near the reported solution, estimate the average angle between (a) the Alg. 2 gradient with K = 100 averaged over many particle-filter runs and (b) a large-K gradient computed with K = 10^4 particles or by finite differences of the log-likelihood estimate. If the K = 100 gradient direction deviates by more than about 30 degrees, or if the larger-K runs change F1 by more than sampling noise, then the K = 100 results do not demonstrate that Eq. (15) minimizes ℓ_R, and the paper needs an explicit bias analysis or a different optimizer.","verdict_should_be":"UNCHANGED","load_bearing_attack":"The paper's central claim—that noisy partial observations recover the polynomial coefficients C and the graph A = 1(abs(C)D^T)—requires that the iteration in Eq. (15) actually minimizes the penalized negative log-likelihood ℓ_R. Two conditions are needed: (i) the gradient estimate from Alg. 2 should be an unbiased, or at least reliably informative, estimator of ∇ℓ, and (ii) stochastic proximal gradient should converge for this nonconvex objective. Neither is established. Reference [13] establishes differentiability of the forward pass, not unbiasedness or variance bounds on the gradient. In Alg. 2, lines 6–7 use a stop-gradient operator on the resampling weights, so the backward pass drops the dependence of the discrete resampling indices on θ; the backpropagated quantity is therefore the gradient of a surrogate, not the gradient of the forward Monte Carlo likelihood estimator. For K = 100 particles in the 20-dimensional chaotic Lorenz 96 system of Sec. V-C, the likelihood estimator is high-variance and the gradient bias can be substantial. Eq. (15) then performs proximal-gradient steps on a stochastic, biased, nonconvex objective; the cited convergence results [43], [44] do not apply as-is. The empirical tables report averages over 150 runs, but the paper does not provide code, the inner-loop step count S, or the tuned λ values, so the reader cannot independently verify that the optimizer reached a good solution rather than a favorable local minimum.","agreement_with_reader":"agree"},"referee_report":{"model":"deepseek-v4-flash","summary":"The paper proposes GraphGrad, a method that approximates the transition density of a general state-space model by a polynomial in the latent state, with a coefficient matrix C and a fixed degree matrix D. The coefficients are estimated by minimizing a penalized negative log-likelihood with an L1 penalty, using gradients obtained from a stop-gradient differentiable particle filter and proximal soft-thresholding updates; a batched variant (B-GraphGrad) is introduced to mitigate likelihood degeneracy. The recovered coefficients define an adjacency graph A = 1_{≠0}(abs(C)D^T), which the authors use to infer directed interactions among state dimensions. The paper reports experiments on Lorenz 63, Lorenz 96, and Kuramoto systems, showing accurate sparse recovery and lower RMSE than a dense polynomial MLE.","tokens_in":25197,"tokens_out":7805,"duration_ms":73639,"significance":"If the claimed performance is robust, GraphGrad fills a useful gap: sparse identification of latent polynomial dynamics under partial and noisy observations. The paper has clear strengths: it uses standard external benchmarks (Lorenz 63, Lorenz 96, Kuramoto), tunes the regularization weight λ on a separate synthetic validation system rather than on the test data, demonstrates scalability to a 20-dimensional Lorenz 96 system with 4,620 parameters, and includes a model-mismatch experiment with an oracle TrueMLE baseline. The reported recovery metrics, including F1 = 1.00 at T = 100 in Table I, are strong if reproducible. The computational cost discussion in Section IV-B is also useful. The main weakness is that the theoretical framing of the optimization is not supported by the provided evidence, and several hyperparameters needed to reproduce the experiments are missing.","major_comments":[{"comment":"The claim that Eq. (15) is a stochastic proximal-gradient method whose convergence is covered by [43], [44] is not supported as stated. The gradient in line 5 of Alg. 4 is obtained from the stop-gradient DPF of Alg. 2, where lines 6-7 stop gradients through the resampling weights. Reference [13] establishes differentiability of the forward pass, but neither unbiasedness of the resulting gradient estimator with respect to ∇ℓ nor a variance bound is shown. For the 20-dimensional Lorenz 96 setting with K = 100, the likelihood estimator is high-variance, and the optimization may be minimizing a surrogate objective rather than ℓ_R in Eq. (12). In addition, Eq. (14) uses the Novograd adaptive update, whereas the convergence results cited for stochastic proximal gradient concern plain proximal-gradient iterations. Please provide a bias/variance analysis or a K-sensitivity study of the recovered C and graph, and qualify the convergence statement accordingly.","section":"III-D2/III-E, Eqs. (14)-(15), Alg. 2"},{"comment":"The central empirical claims are not fully reproducible from the manuscript because key training hyperparameters are not reported. Section V-A fixes B = ceil(T/10), η = 10^-3, and K = 100, but it does not give the number S of inner optimization steps per batch, and the final tuned λ values are never reported. Section V-D does not state the polynomial degree d used for the Kuramoto experiment. The paper also does not provide code or data repositories. Please report S and the tuned λ for every table and figure, state d for the Kuramoto experiment, and make code and data available or provide a complete hyperparameter table.","section":"V-A, V-B1, V-C, V-D"},{"comment":"The definition and reporting of RMSE need clarification. Section V-A defines numerically zero as |C_ij| < 10^-6 and says RMSE measures entry recovery, while Section V-B1 states that RMSE is computed only on terms recovered as non-zero; these are different quantities. Table IX reports pMLE RMSE 0.32 at T = 25, which is inconsistent with the surrounding values (pMLE is otherwise 1.3-2.5 in that table) and with the text claiming that pMLE performs poorly. Please verify this entry and state precisely how RMSE is computed in each table.","section":"V-A, Table IX"}],"minor_comments":[{"comment":"The abstract contains an incomplete phrase, '...a rich and efficient d of a general SSM', which should be corrected.","section":"Abstract"},{"comment":"The weight update writes 'w(i)_t' where the index should presumably be the particle index k; please correct the notation.","section":"Alg. 1 and Alg. 2, line 10"},{"comment":"The tables report only point averages over 150 runs for the classification metrics; please add standard deviations or intervals, as is done for RMSE in the figures.","section":"Section V-A"},{"comment":"Reference [10] is spelled 'SINdY' in the text; the correct name is SINDy.","section":"Section IV-D and Refs"},{"comment":"For the d = 3 Lorenz 96 results, pMLE at T = 25 is reported as 0.32 yet the same method is reported as 2.5 at T = 50; please check whether this is a typo or an artifact of the RMSE definition.","section":"Table IX"},{"comment":"The caption says red links correspond to links found with an incorrect monomial, but the figure legend in the right plot is not explicit; please make the blue/red distinction clear in the printed legend.","section":"Figure 5"}],"recommendation":"major_revision","confidential_remarks":"The empirical study is substantial and the methodology is a sensible contribution to sparse estimation in nonlinear SSMs. My main reservation is that the convergence and gradient properties are asserted more strongly than they are established, and the missing hyperparameters and artifacts make the headline F1 = 1.00 results hard to verify. I would support acceptance after the authors either add the missing analysis and details or explicitly weaken the theoretical claims."},"author_rebuttal":null,"desk_editor":{"model":"deepseek-v4-flash","letter":"Colleague,\n\nRecommended for review. This is a genuine improvement over the authors' ICASSP predecessor, and the experiments stand up to inspection, but the optimization story is softer than the text implies and the reproducibility details are incomplete.\n\nWhat's new: they write the transition as a fixed degree matrix D and a learnable coefficient matrix C, so the graph structure comes from abs(C)D^T. They train C through a differentiable particle filter with stop-gradient resampling, then apply an L1 proximal step. The proximal update (rather than subgradient) and the telescoping batches are the actual additions over [16]. The batching argument is sensible: fit a coarse model on short prefixes and refine, which mitigates likelihood degeneracy in long series.\n\nWhat's good: the simulations are clean. Lorenz 63 at d=2 recovers structure to F1=1.00 by T=100, and Lorenz 96 (20 dimensions, 4620 parameters) reaches the same by T=200. The Kuramoto model-mismatch experiment is a useful stress test; nRMSE ~1.2 versus the oracle TrueMLE is a respectable result. The tuning of lambda on a separate synthetic system is appropriate. The runtime scaling discussion and the parallelisation remark are honest.\n\nSoft spots: the convergence statement leans on [43],[44], which are convex results, while the likelihood is nonconvex and the DPF gradient is biased. The stop-gradient in Alg. 2 drops the dependence of resampling indices on C, so Eq. (15) is not, strictly, a proximal-gradient iteration on the penalized log-likelihood. This is a real gap, but the empirical recovery is consistent, so I would not call it fatal. The bigger problem is reproducibility: the paper does not report the inner-loop step count S, the tuned lambda values, or provide code. Without those, the tables cannot be independently checked. The baselines are also narrow—only the dense pMLE and the TrueMLE oracle—though SINDy and relatives need direct state observations, so the comparison set is defensible.\n\nWhom is this for? Signal-processing and statistics people working on parameter estimation in NLSSMs, and practitioners in ecology or epidemiology who want interpretable connectivity from noisy partial observations. It deserves a serious referee. Send it, but with a request for code, S, lambda values, a K-sensitivity check, and a softer statement on convergence guarantees.","headline":"A solid engineering contribution that combines a differentiable particle filter with proximal L1 for sparse recovery of polynomial SSM dynamics; referee it, but ask for reproducibility details and a more careful statement of the optimization guarantees.","tokens_in":25784,"tokens_out":4565,"would_cite":true,"duration_ms":37960,"reading_group":"yes","serious_thinker":"yes","would_accept_peer_review":true},"rs_alignment":null,"lean_confirmation":null,"pith_extraction":{"msc":["62M09","62M05","65K10"],"pacs":[],"model":"deepseek-v4-flash","headline":"From noisy, partial observations, a sparse polynomial model can recover the interaction graph of a nonlinear dynamical system.","keywords":["state-space models","differentiable particle filter","sparse polynomial approximation","nonlinear dynamics","network inference","L1 regularization","proximal gradient","Lorenz systems"],"falsifier":"Simulate a sparse polynomial SSM with known coefficients, run B-GraphGrad on very long observation series with a correctly specified degree, and check whether the recovered adjacency matrix A converges to the true support; the paper's claim predicts F1 approaching 1, and any plateau below 1 at large T refutes it.","tokens_in":24663,"feed_emoji":"🔗","tokens_out":9353,"duration_ms":74219,"temperature":0.7,"pith_summary":"GraphGrad tackles a standing problem in time-series analysis: when a system's hidden state is only seen through noisy observations, how can one learn which components actually influence one another? The paper's answer is to approximate the unknown transition function by a polynomial of fixed maximum degree, with coefficients stored in a matrix C, and to read the zero pattern of C as a directed interaction graph. It shows that a differentiable particle filter can supply gradients of a Monte Carlo likelihood estimate, and that a proximal step applying L1 sparsity turns those gradients into a practical estimator for C. The result, if correct, is a fully automatic, unsupervised route from raw observation series to both the structure and the parameter values of a broad class of nonlinear state-space models.","feed_headline":"Sparse polynomial fits expose interaction graphs","feed_subtitle":"A differentiable particle filter turns noisy observations into sparse transition coefficients and a directed interaction graph.","key_machinery":"The load-bearing object is the pair (C, D): a learned real coefficient matrix C ∈ $R^{{Nx×M}}$ and a fixed integer degree matrix D ∈ $N_0^{{Nx×M}}$ whose columns list the exponents of all monomials of degree at most d in N_x variables. The transition mean is the polynomial f(x, C; D) = Σ_{j=1}^M C_{·,j} ∏_{i=1}^{N_x} $x_i^{{D_{i,j}}$}. Gradient information flows to C through a stop-gradient differentiable particle filter (Algorithm 2), which leaves the forward particle filter unchanged but lets automatic differentiation pass through the importance weights; the L1 penalty is then applied not by subgradient descent but by the soft-thresholding proximity operator T_{ηλ}(x) = max(|x|−ηλ, 0) sgn(x), giving the stochastic proximal-gradient update of Equation (15). A telescoping batch schedule (B-GraphGrad) avoids likelihood degeneracy on long observation series.","core_discovery":"The central claim is that, for a nonlinear state-space model whose transition density is written as a Gaussian around a polynomial mean f(x, C; D) = Σ_j C_{·,j} ∏_i $x_i^{{D_{i,j}}$}, the coefficient matrix C can be estimated from the observations alone by minimizing a penalized negative log-likelihood, with the L1 penalty applied through a soft-thresholding proximal update. Because the degree matrix D is fixed and known, every nonzero coefficient names a specific monomial coupling one state to another, so the matrix A = 1_{≠0}(|C|D^T) defines a directed graph of state interactions. The paper reports that this scheme recovers the connectivity of the Lorenz 63 and Lorenz 96 systems with near-perfect accuracy, and that on the non-polynomial Kuramoto oscillator it comes within about 20% of the state-recovery error of an oracle that knows the true model form.","pith_inferences":["The reliability of the recovered graph should fall as observation noise grows and as the polynomial degree increases; a natural extension is to report support-recovery F1 as a function of signal-to-noise ratio, which the paper only partially tabulates.","Because the estimator is a function-library regression, swapping monomials for trigonometric or exponential basis terms would let the same pipeline handle oscillatory systems; the paper notes this possibility but does not test it.","The differentiable particle filter gradient is not proven unbiased, so practitioners should treat recovered coefficients as point estimates and could diagnose instability by comparing multiple random restarts.","For very high-dimensional systems the number of monomials M grows combinatorially; the paper's complexity discussion covers parallelization but not whether the statistical estimator remains practical at that scale."],"forward_implications":["For any system whose transition function is polynomial, the method yields both the coefficient values and the zero-pattern graph from observations only, without any state labels.","The recovered graph A = 1_{≠0}(|C|D^T) provides a directed, Granger-style notion of which state dimensions drive which.","Under model mismatch, the polynomial surrogate still supports state estimation; on the Kuramoto oscillator the paper reports state-recovery error only about 20% above an oracle that knows the true model.","The telescoping batch schedule removes the likelihood-degeneracy obstacle, so long observation series can be used without vanishing gradients.","The approach extends to 20-dimensional states and to over-specified polynomial degrees with only modest loss in support recovery."],"supporting_citations":[{"why":"Supplies the stop-gradient differentiable particle filter that makes gradient-based likelihood optimization possible.","marker":"[13]"},{"why":"Provides the Monte Carlo likelihood estimator and particle-filter background on which the loss function rests.","marker":"[17]"},{"why":"Justifies the proximal soft-thresholding update for the L1 penalty.","marker":"[42]"},{"why":"Cited for the convergence theory of stochastic proximal gradient methods underpinning the optimization.","marker":"[43]"},{"why":"The SINDy sparse-identification method that motivates the function-library regression interpretation.","marker":"[10]"},{"why":"The Lorenz 63 system serves as the main exactly polynomial test case.","marker":"[14]"},{"why":"The Lorenz 96 system provides the high-dimensional chaotic test case.","marker":"[15]"},{"why":"The Kuramoto oscillator is the model-mismatch test system that cannot be exactly represented.","marker":"[49]"}],"fun_headline_variants":["Sparse polynomials reveal hidden system graphs","GraphGrad: Sparse fits unmask state interactions","Polynomial sparsity uncovers dynamic coupling","Particle filter learns interaction networks","Sparse SSM fitting exposes directed graphs"],"cache_read_input_tokens":3200,"weakest_assumption_plain":"The whole pipeline works only if the noisy gradient returned by the differentiable particle filter, computed from a modest number of particles, is a reliable enough descent direction in a nonconvex landscape; the paper does not prove this, and a biased or trapped gradient could yield a wrong graph even with unlimited observations.","fun_headline_variants_meta":{"raw":{"variants":["Sparse polynomials reveal hidden system graphs","GraphGrad: Sparse fits unmask state interactions","Polynomial sparsity uncovers dynamic coupling","Particle filter learns interaction networks","Sparse SSM fitting exposes directed graphs"]},"model":"deepseek-v4-flash","effort":"low","cost_usd":0.000169,"raw_usage":{"total_tokens":1288,"prompt_tokens":993,"completion_tokens":295,"prompt_tokens_details":{"cached_tokens":384},"prompt_cache_hit_tokens":384,"prompt_cache_miss_tokens":609,"completion_tokens_details":{"reasoning_tokens":229}},"tokens_in":609,"tokens_out":295,"duration_ms":3395,"temperature":1.0,"reasoning_tokens":229,"cache_read_input_tokens":384,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-08-12T14:05:12.572015+00:00","model_set":{"reader":"deepseek-v4-flash"},"falsifier":"Simulate a sparse polynomial SSM with known coefficients, run B-GraphGrad on very long observation series with a correctly specified degree, and check whether the recovered adjacency matrix A converges to the true support; the paper's claim predicts F1 approaching 1, and any plateau below 1 at large T refutes it.","supporting_citations":[{"cited_title":"S¨arkk¨a, Bayesian Filtering and Smoothing","cited_arxiv_id":null,"evidence_quote":"Provides the Monte Carlo likelihood estimator and particle-filter background on which the loss function rests."},{"cited_title":"Proximal splitting methods in signal processing,","cited_arxiv_id":null,"evidence_quote":"Justifies the proximal soft-thresholding update for the L1 penalty."},{"cited_title":"Convergence of stochastic proximal gradient algorithm,","cited_arxiv_id":null,"evidence_quote":"Cited for the convergence theory of stochastic proximal gradient methods underpinning the optimization."},{"cited_title":"Deterministic nonperiodic flow,","cited_arxiv_id":null,"evidence_quote":"The Lorenz 63 system serves as the main exactly polynomial test case."},{"cited_title":"Predictability: A problem partly solved,","cited_arxiv_id":null,"evidence_quote":"The Lorenz 96 system provides the high-dimensional chaotic test case."},{"cited_title":"Kuramoto,Chemical Oscillations, Waves, and Turbulence","cited_arxiv_id":null,"evidence_quote":"The Kuramoto oscillator is the model-mismatch test system that cannot be exactly represented."}],"review_version":1}