{"id":"b95c24c6-4dfe-4af7-9b7a-159e077162c6","arxiv_id":"2506.15906","paper_version":1,"verdict":"CONDITIONAL","confidence":"MODERATE","novelty_score":5.0,"correctness_risk":"medium","formal_verification":"none","parameter_count":4,"one_line_summary":"LoGoS-GPO combines nearest-neighbor spatial kernels, inducing-point sparsity, Kronecker structure, and a wavelet neural operator mean to make Gaussian Process Operators scalable while retaining uncertainty quantification.","lead":"This paper proposes LoGoS-GPO, a faster Gaussian process model for learning physical systems from data, using local sparsity and structured kernels to cut computation time. A generalist might care because it aims to combine the accuracy of neural operator models with reliable uncertainty estimates for large-scale simulation.","discovery_kind":"new_method","skeptic_critique":{"model":"deepseek-v4-flash","headline":"The KNN-truncated spatial kernel in Eq. (13) is not guaranteed positive definite, so the inverse in Eq. (20) and the claimed O(dK) complexity do not follow from the stated construction.","rationale":"The paper's central claim is that the combination of Kronecker factorization, KNN spatial sparsification, and sparse variational inference reduces GPO cost from O(N^3 d^3) to O(M^3 + M^2 B + dK^2) while maintaining accuracy. I read Section 3.1 as requiring that K_N^x in Eq. (13) be a valid covariance matrix whose inverse can be exploited through Eq. (20). The reader's weakest-assumption analysis identified exactly this point, and my closer reading sharpens it: KNN truncation is not PD-preserving, so the matrix can be indefinite in exactly the fine-resolution regimes the paper targets. Even in PD cases, the inverse of a sparse matrix is generally dense, so the claimed inversion cost needs an additional argument about sparse factorization or approximate inversion. The empirical results in Section 4 are internally consistent and the accuracy numbers are plausible, but they do not resolve the mathematical issue because K and the kernel hyperparameters are not reported, and the code is withheld until acceptance. The missing ablations of the enhanced WNO mean and the non-zero mean are secondary; the PD/inverse concern is load-bearing. The construction is repairable in principle, for example with jitter or a Vecchia-style formulation, so I would not escalate the verdict to rejection. The conditional verdict already given by the reader remains appropriate, so no change is needed.","tokens_in":19174,"tokens_out":6000,"duration_ms":67295,"concrete_test":"Reproduce Algorithm 1 on a uniform periodic 1D grid with d = 1024, using the RBF kernel k(x, x') = exp(-||x - x'||^2 / (2 l^2)) with l equal to the grid spacing and, if available, the actual hyperparameters used in the Burgers experiment. Build the KNN-truncated matrix P by retaining K = 5 nearest neighbors per row (including the diagonal), symmetrize it as (P + P^T) / 2, and compute its smallest eigenvalue. Repeat for K = 10, 20, 50 and for l / Delta_x in {0.5, 1, 2}. If any eigenvalue is negative, Eq. (13) is not a valid covariance matrix and Eq. (20) cannot be used as written; the authors should then add jitter, adopt a PD-preserving local approximation, or revise the complexity claim. As a secondary output, for the cases that are positive definite, compute the exact inverse and count nonzeros versus dK to test the O(dK) inverse-density assumption.","verdict_should_be":"UNCHANGED","load_bearing_attack":"The central scalability claim in Section 3.1 rests on treating the KNN-truncated matrix K_N^x from Eq. (13) as a covariance matrix whose inverse can be used in Eq. (20) at O(dK) or O(dK^2) cost. This is not justified by the construction. KNN truncation does not preserve positive definiteness: on a uniform periodic 1D grid with an RBF kernel and fine spacing, the resulting banded matrix has off-diagonal correlations that can approach 1, giving negative eigenvalues. For a tridiagonal section with correlation rho, the smallest eigenvalue is 1 - 2*rho (up to boundary terms), which is negative whenever adjacent-point correlation exceeds 0.5. If K_N^x is indefinite, the matrix inverse in Eq. (20) is not a valid covariance inverse, the Kronecker identity is undefined, and the ELBO in Eq. (21) is not a valid lower bound. Second, even when K_N^x is positive definite, the inverse of a sparse matrix is generally dense; the claimed O(dK) inversion cost requires either sparse Cholesky factorization with bounded fill-in or an explicit sparse-inverse approximation, and neither is stated or derived. Additionally, the KNN neighborhood relation is not symmetric in general, so the raw matrix from Eq. (13) need not even be symmetric, although a covariance matrix must be. The construction could be repaired with symmetrization, a diagonal jitter, or a Vecchia-style sequential conditioning that guarantees a valid sparse precision matrix, but none of these is present. Without this, the theoretical complexity advantage and the Kronecker-SVGP inference steps do not follow from the stated method, and the empirical speedups cannot be attributed to the claimed mechanism.","agreement_with_reader":"agree"},"referee_report":{"model":"deepseek-v4-flash","summary":"This paper proposes LoGoS-GPO, a sparse variational Gaussian process operator that combines a Kronecker-factorized covariance over (input-feature, spatial-grid) pairs, K-nearest-neighbor truncation of the spatial covariance, and a neural-operator mean and latent embedding. The central scaling claim in Section 3.1 is that this reduces training cost to O(M^3 + M^2B + dK^2) from O(N^3d^3) for exact inference. The method is evaluated on Burgers, wave advection, Darcy, and Navier-Stokes benchmarks, reporting relative L2 errors of 0.43%-2.01%, which are better than GPO and SVGP, together with scaling experiments in grid size and sample count. An enhanced wavelet neural operator is introduced to provide the mean function and the latent kernel embedding.","tokens_in":19526,"tokens_out":12431,"duration_ms":126720,"significance":"If the scalability claim were established, LoGoS-GPO would be a practically useful step toward uncertainty quantification in operator learning, and the empirical comparison against SVGP and GPO on standard PDE benchmarks is a genuine strength. The paper makes a concrete, falsifiable algorithmic proposal (KNN sparsification plus Kronecker factorization plus sparse variational inference) and reports external benchmarks rather than constructing the evaluation from its own prior outputs. However, the advertised complexity reduction rests on unproved spectral and algorithmic properties of the KNN-truncated spatial matrix, and the current derivation does not support the central claim as written. The empirical results are plausible but cannot by themselves validate the complexity model.","major_comments":[{"comment":"The KNN-truncated matrix K_N^x defined in Eq. (13) is not guaranteed to be symmetric, because the K-nearest-neighbor relation is not symmetric, and it is not guaranteed to be positive definite. Even for a symmetrized banded version, truncating an RBF kernel on a fine uniform grid gives a tridiagonal section whose smallest eigenvalue is approximately 1 - 2ρ, which is negative whenever the adjacent-point correlation ρ exceeds 0.5. Since Eq. (20) uses (K_N^x)^{-1} in a Kronecker inverse and Eqs. (21)-(23) treat K_N^x as a covariance, the variational objective and predictive variance are not valid as stated. The construction needs to be repaired, for example by symmetrization with a diagonal jitter or by replacing the KNN truncation with a Vecchia-style sparse precision that is guaranteed positive definite.","section":"§3.1, Eq. (13) and Eq. (20)"},{"comment":"Algorithm 1 computes all pairwise distances for each grid point, so constructing K_N^x costs O(d^2) kernel evaluations and O(d^2 log d) sort time, which contradicts the stated O(d log N + dK) construction and inversion cost. Furthermore, the inverse of a sparse matrix is generally dense, and no argument is given that K_N^x admits a bounded-fill-in Cholesky factorization or that its inverse can be applied in O(dK). The overall complexity claim O(M^3 + M^2B + dK^2) is therefore not established; please specify the exact linear algebra (bandwidth, fill-in, or iterative solver) and derive the cost from the actual algorithm.","section":"§3.1, Algorithm 1 and complexity statement"},{"comment":"The derivation of the enhanced wavelet convolution is not self-contained. Equation (28) gives a convolution in the wavelet translation variable t, but Eq. (29) replaces it with F^{-1}[R_φ · F(Wv_j)(s,ω)] without specifying the domain of the Fourier transform, the relation between R_φ and K_φ, or how the discrete wavelet transform and FFT implementations compose. Please provide a precise statement of the transforms and variables, or give the actual computational formula used in the experiments; without this, the claimed improvement over the original WNO cannot be verified.","section":"§3.2, Eqs. (27)-(31)"},{"comment":"Equation (12) asserts the separable Kronecker decomposition K = K_a ⊗ K_x, but this does not follow from the latent-space kernel in Eq. (9) unless Φ and k_ax are assumed to factorize over the input-feature and spatial coordinates. The paper should state this explicitly as a modeling assumption and indicate what class of kernels is covered; as written, the decomposition is introduced without derivation or justification.","section":"§3.1, Eq. (12)"}],"minor_comments":[{"comment":"Equation (25) uses g*(ω) for the Fourier transform and integrates C_g over Ω with respect to ω; please standardize the notation (e.g., \\hat g or \\mathcal{F}g) and state the integration domain as the frequency axis.","section":"§3.2, Eq. (25)"},{"comment":"The periodic boundary condition u(x - π) = u(x + π) on x ∈ (0,1) is inconsistent with the stated domain; it should presumably read u(0,t) = u(1,t).","section":"§4.2, Eq. (34)"},{"comment":"The learning rates are reported in an ambiguous product notation such as '58×10^-3' and '32×10^-2'; please use standard scientific notation.","section":"Table 2"},{"comment":"Algorithm 2 reuses the symbol B for both the mini-batch subset and the mini-batch size; rename one of them to avoid confusion.","section":"Algorithm 2"},{"comment":"The text says Fig. 8 investigates sample-size variation, but Fig. 8 shows model predictions; the reference should probably be to Fig. 9b.","section":"§4.4"},{"comment":"The relative L2 error used in Table 3 is not defined; please provide the formula used to compute it.","section":"Table 3"},{"comment":"The notation for kernel hyperparameters is inconsistent: Section 2.2 defines θk = {θ_{k1}, θ_{k2}}, while Section 3.1 uses θk = {θ_{k3}, θa, θx}; please unify the notation.","section":"§3.1 vs. §2.2"},{"comment":"The paper claims the enhanced WNO as a contribution, but no experiment compares LoGoS-GPO using the enhanced WNO with LoGoS-GPO using the original WNO; without such an ablation, the incremental benefit of the enhanced WNO is not demonstrated.","section":"§4, enhanced WNO"}],"recommendation":"major_revision","confidential_remarks":"The central scalability claim needs a rigorous repair before the paper can be accepted. The KNN sparsification as written is not a valid covariance approximation; replacing it with a Vecchia-style construction, symmetrization plus jitter, or a bounded-fill-in factorization is likely feasible but would require re-running the experiments and re-deriving the complexity. The enhanced-WNO contribution is not isolated experimentally, and the manuscript leans heavily on the authors' prior GPO and WNO work; a clearer statement of incremental novelty would help the editors assess fit. The empirical comparison against GPO and SVGP on standard PDE benchmarks is a strength and should be preserved."},"author_rebuttal":null,"desk_editor":{"model":"deepseek-v4-flash","letter":"First thing you should know: this paper is worth reading before you cite it for scalable Gaussian Process Operators. The combination of KNN spatial sparsification, SVGP in feature space, Kronecker factorization, and a neural-operator mean is genuinely new, and the empirical results on four PDE benchmarks are strong enough that the method will probably work if the theory is fixed. But the load-bearing complexity claim does not hold as written.\n\nThe new integration is real. Eq. (12) factorizes the covariance into feature and spatial Kronecker factors, Eq. (13) sparsifies the spatial factor by K-nearest neighbors, and Eq. (21) uses SVGP over inducing points. The enhanced WNO in Section 3.2 also identifies a genuine limitation of the original WNO (wavelet-domain elementwise multiplication is not exact convolution) and proposes a Fourier-filtered wavelet convolution. The experiments show 0.43–2.01% relative L2 error, beating GPO and SVGP, with wall-clock and memory advantages that look monotone in grid size and sample count. That is a useful contribution to the UQ-for-PDEs discussion.\n\nThe soft spots are in the theory, and one is serious. The KNN-truncated matrix from Algorithm 1 is not symmetric as constructed, because membership in K nearest neighbors is not reciprocal. Even after symmetrization, sparsification can destroy positive definiteness: on a uniform grid with close neighbors, the off-diagonal correlations can approach 1, and a tridiagonal section can have negative eigenvalues. The paper uses [K_x]^{-1} in Eq. (20) and claims an O(dK) inversion cost, but the inverse of a sparse matrix is generally dense. To get an O(dK) inverse you need a sparse Cholesky with bounded fill-in or a Vecchia-style sparse precision matrix; neither is stated. Without that, the ELBO in Eq. (21) is not a valid lower bound and the complexity claim collapses. This is fixable with symmetrization plus jitter or a proper sparse precision construction, but the paper needs to do it.\n\nTwo smaller issues. The enhanced WNO derivation in Eqs. (27)–(30) mixes wavelet and Fourier variables; the final form may be correct, but the derivation is too terse to verify. And there is no ablation comparing enhanced WNO to the original WNO, so we cannot attribute the accuracy gains to that change. The paper also withholds code until acceptance and does not report K or M or any calibration metric, so \"calibrated uncertainty\" is asserted, not shown.\n\nWho is this for? Anyone working on scalable GP-based operator learning. It deserves a serious referee, but the referee should first ask for a valid sparse-inverse argument and the ablations. I would not cite it in its current form.","headline":"A promising blend of sparsification tricks for GPOs, but the central complexity claim rests on an unproven and likely invalid assumption about the inverse of the KNN-truncated kernel.","tokens_in":20093,"tokens_out":3709,"would_cite":false,"duration_ms":36080,"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":"This paper proposes LoGoS-GPO, a Gaussian Process Operator that combines K-nearest-neighbor spatial kernels, Kronecker factorization, and sparse variational inference to lower operator-learning inference cost from $O(N^3 d^3)$ to $O(M^3 +…","keywords":["Gaussian process operator","operator learning","Kronecker factorization","K-nearest-neighbor kernel","sparse variational inference","wavelet neural operator","uncertainty quantification","parametric PDEs"],"falsifier":"Compute the smallest eigenvalue of the truncated matrix $K_x$ from Eq. (13) on a uniform 1D grid with an RBF kernel and small $K$; if any eigenvalue is negative for typical lengthscales, the matrix is not a valid covariance and the posterior variance in Eq. (23b) can go negative. Separately, benchmark the inversion step's wall-clock scaling in $d$ at fixed $K$: growth like $d^2$ or $d^3$ rather than $O(dK)$ would falsify the headline complexity claim.","tokens_in":18980,"feed_emoji":"📐","tokens_out":9052,"duration_ms":75807,"temperature":0.7,"pith_summary":"This paper argues that a Gaussian Process Operator—a probabilistic model for learning maps between function spaces—can be made scalable without sacrificing accuracy by combining three structural shortcuts: a Kronecker product splitting the kernel into parameter-space and spatial parts, a K-nearest-neighbor truncation of the spatial covariance, and sparse variational inference over inducing points in the parameter space. On Burgers, wave advection, Darcy flow, and Navier–Stokes benchmarks, the method reports relative L2 errors between 0.43% and 2.01%, lower than the original GPO and a standard sparse GP baseline. If the claim holds, it removes the cubic-complexity barrier that has kept probabilistic neural operators off high-resolution PDE problems, while preserving predictive uncertainty estimates. The paper also introduces an enhanced wavelet neural operator that supplies a non-zero prior mean and a latent-space kernel, which is what the authors credit for compensating the accuracy lost to local sparsification.","feed_headline":"Local sparse Gaussian process operator cuts PDE cost to near-linear","feed_subtitle":"GPOs combine local spatial kernels, Kronecker structure, and sparse inference to beat cubic scaling on nonlinear PDEs.","key_machinery":"The load-bearing object is the doubly sparse covariance $K = K_a(\\Phi(A),\\Phi(A')) \\otimes K_x(x,x')$: a Kronecker product of a parameter-space kernel evaluated on neural-operator embeddings and a spatial kernel $K_x$ truncated so each row retains only its $K$ nearest neighbors. Kronecker structure lets matrix products and inverses act on the two factors separately, and sparse variational inference with $M$ inducing points replaces the full $Nd$-dimensional posterior with a low-rank one, giving the claimed $O(M^3 + M^2B + dK^2)$ cost. The non-zero mean $m(a,x)$, produced by an enhanced wavelet neural operator, is not decorative: it carries the global, operator-aware information that the localized spatial kernel discards, and the same neural operator defines the latent space $\\Phi(\\cdot)$ in which the parameter-space kernel is evaluated. The paper's formal route from function-valued GP to this finite model is the probabilistic currying equivalence $G(a)(x) = f(a,x)$.","core_discovery":"The central claim is that a doubly sparse Gaussian Process Operator can match or beat exact GPO accuracy at a fraction of the cost. The proposed model writes the covariance as $K = K_a(\\Phi(A),\\Phi(A')) \\otimes K_x(x,x')$, keeps only the $K$ nearest neighbors of each spatial grid point in $K_x$, and treats the parameter-space factor $K_a$ with sparse variational inference using $M$ inducing points. The paper states that this combination reduces inference cost from $O(N^3 d^3)$ to $O(M^3 + M^2B + dK^2)$, where $N$ is the number of training samples, $d$ the discretization size, and $B$ the mini-batch size. On four nonlinear PDE benchmarks—1D Burgers, 1D wave advection, 2D Darcy flow on a notched triangular domain, and 2D Navier–Stokes—the reported relative L2 errors are 0.86%, 0.43%, 1.38%, and 2.01%, respectively, all below the corresponding GPO and SVGP errors. A separate claim is that the enhanced wavelet neural operator, which replaces the pointwise-multiplication approximation of wavelet convolution with a true convolution integral in the wavelet domain, provides the expressive mean and latent embedding that let the sparsified kernel stay accurate.","pith_inferences":["The step most worth stress-testing is the claimed $O(dK)$ inversion of the KNN-truncated spatial matrix: a sparse matrix with $K$ nonzeros per row is not guaranteed to be positive definite, and its inverse is generally dense, so the stated cost depends on an extra structural assumption the paper leaves implicit.","The accuracy comparisons are against GPO and SVGP; a sharper test of whether sparsification costs accuracy would compare against deterministic wavelet or Fourier neural operators at the same resolution and sample budget.","A natural extension is to make the neighbor set adaptive to kernel lengthscales or learned, which could recover long-range spatial interactions while retaining the sparse layout.","One could probe the role of the operator-aware mean directly by ablating it: if accuracy collapses when the mean is set to zero, the expressive prior mean is doing the work that the sparse kernel alone cannot."],"forward_implications":["If the complexity claim holds, GPOs become trainable at grid resolutions and sample sizes that were previously out of reach, with per-epoch time and memory growing roughly linearly in the grid size $d$ rather than cubically.","The four benchmarks suggest that probabilistic operator learning does not have to trade uncertainty quantification for scalability: predictive intervals remain narrow on unseen inputs from the same distribution, including near the discontinuity in the wave advection case.","The non-zero neural-operator mean is a transferable recipe: any discretization-invariant operator could supply the mean and the latent embedding, turning a deterministic neural operator into a scalable Gaussian process with uncertainty estimates.","Replacing the pointwise-multiplication approximation in wavelet convolution with a genuine convolution integral in the wavelet domain is presented as an architectural fix that preserves multiscale locality while regaining global interactions.","If resolution independence transfers as reported, the model offers a route to training on coarse grids and predicting on finer grids with calibrated uncertainty, which is the practical setting for parametric PDE surrogates."],"supporting_citations":[{"why":"Defines the original Gaussian Process Operator and provides the baseline whose exact inference cost the paper aims to reduce.","marker":"[21]"},{"why":"Supplies the sparse variational inference scheme over inducing points used for the parameter-space kernel.","marker":"[22]"},{"why":"Supplies the nearest-neighbor local kernel approximation strategy that the paper adapts to the spatial covariance.","marker":"[24]"},{"why":"Contributes the Kronecker-structured covariance machinery that motivates the separable kernel factorization.","marker":"[27]"},{"why":"Introduces the wavelet neural operator architecture that the paper enhances and uses for the mean and latent embedding.","marker":"[4]"},{"why":"Provides the probabilistic currying equivalence that lets the function-valued GP be treated as a finite vector-valued GP.","marker":"[33]"},{"why":"Supplies the Burgers' equation dataset and experimental setup used in the first case study.","marker":"[3]"}],"fun_headline_variants":["Sparse kernel trick tames cubic GPO cost on nonlinear PDEs","Nearest-neighbor kernels make Gaussian process operators scale","Local kernels and Kronecker factors beat cubic GPO scaling","Scalable GPO: local kernels, sparse inference, near-linear cost","Doubly sparse GPO matches full accuracy at fraction of cost"],"cache_read_input_tokens":3200,"weakest_assumption_plain":"The load-bearing premise is that the K-nearest-neighbor-truncated spatial covariance matrix $K_x$ is still a valid covariance whose inverse can be computed in $O(dK)$ time; a sparse matrix with $K$ nonzeros per row is not generally positive definite, and its inverse is generally dense, so Eq. (20)'s use of $[K_x]^{-1}$ and the quoted $O(dK)$ inversion cost are not automatic.","fun_headline_variants_meta":{"raw":{"variants":["Sparse kernel trick tames cubic GPO cost on nonlinear PDEs","Nearest-neighbor kernels make Gaussian process operators scale","Local kernels and Kronecker factors beat cubic GPO scaling","Scalable GPO: local kernels, sparse inference, near-linear cost","Doubly sparse GPO matches full accuracy at fraction of cost"]},"model":"deepseek-v4-flash","effort":"low","cost_usd":0.000307,"raw_usage":{"total_tokens":1809,"prompt_tokens":1050,"completion_tokens":759,"prompt_tokens_details":{"cached_tokens":384},"prompt_cache_hit_tokens":384,"prompt_cache_miss_tokens":666,"completion_tokens_details":{"reasoning_tokens":671}},"tokens_in":666,"tokens_out":759,"duration_ms":6766,"temperature":1.0,"reasoning_tokens":671,"cache_read_input_tokens":384,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-08-15T19:29:49.910670+00:00","model_set":{"reader":"deepseek-v4-flash"},"falsifier":"Compute the smallest eigenvalue of the truncated matrix $K_x$ from Eq. (13) on a uniform 1D grid with an RBF kernel and small $K$; if any eigenvalue is negative for typical lengthscales, the matrix is not a valid covariance and the posterior variance in Eq. (23b) can go negative. Separately, benchmark the inversion step's wall-clock scaling in $d$ at fixed $K$: growth like $d^2$ or $d^3$ rather than $O(dK)$ would falsify the headline complexity claim.","supporting_citations":[{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"Defines the original Gaussian Process Operator and provides the baseline whose exact inference cost the paper aims to reduce."},{"cited_title":"Variational nearest neighbor gaussian process","cited_arxiv_id":null,"evidence_quote":"Supplies the nearest-neighbor local kernel approximation strategy that the paper adapts to the spatial covariance."},{"cited_title":"Linearization turns neural operators into function-valued gaussian processes, 2024","cited_arxiv_id":null,"evidence_quote":"Provides the probabilistic currying equivalence that lets the function-valued GP be treated as a finite vector-valued GP."}],"review_version":1}