{"id":"78e846fb-16fa-4dd0-a380-80c5b2c0e9ae","arxiv_id":"2507.06677","paper_version":1,"verdict":"CONDITIONAL","confidence":"MODERATE","novelty_score":6.0,"correctness_risk":"medium","formal_verification":"none","parameter_count":2,"one_line_summary":"RLRTO-based sampling enforces monotonicity constraints in Gaussian processes with independent samples and large speedups over Gibbs and NUTS samplers.","lead":"This paper proposes using an optimization-based sampler, regularized linear randomize-then-optimize (RLRTO), to enforce monotonicity constraints on Gaussian process surrogates, replacing slower MCMC samplers. It reports comparable predictive accuracy to existing constrained GP methods with substantially higher computational efficiency, demonstrated on synthetic functions and differential-equation surrogates.","discovery_kind":"new_application","skeptic_critique":{"model":"deepseek-v4-flash","headline":"RLRTO produces a censored (projected) distribution with point mass at zero, not the truncated-Gaussian posterior targeted by the baselines; without a test of this implicit posterior, the comparison in §4 is between different models, not samplers.","rationale":"The paper's main contribution is a fast sampling scheme for monotonicity-constrained GPs. The method's validity rests on the claim that RLRTO draws from a constrained posterior. A careful reading of the algorithm shows that the constrained solution is the projection (in the Hessian norm) of a draw from the unconstrained Gaussian posterior onto the non-negative orthant. For a scalar model this is the censored variable max(0, Z), which assigns positive probability to the boundary zero. The truncated-Gaussian prior method's posterior (Eq. 29) has no boundary atom. Hence, RLRTO and the baselines are different generative models. The paper acknowledges the density is implicit, but the abstract's phrasing and the way results are compared (same MSE/CI metrics, claims of 'comparable predictive performance') imply a common inference target. This is misleading: if the RLRTO posterior shrinks derivatives to zero much more aggressively, its CI widths and predictive accuracy are not directly comparable to the truncated methods. A concrete check is to measure the boundary activation rate and compare predictive intervals against the truncated-Gaussian posterior in one experiment; if they differ, the central claim must be rephrased. This does not invalidate the method as a heuristic constrained GP, but it does invalidate the unqualified 'posterior sampling' claim and makes the experimental comparison a comparison of models, not samplers. The reader's conditional verdict is appropriate; no further adjustment is needed.","tokens_in":19411,"tokens_out":15179,"duration_ms":167740,"concrete_test":"Take the 1D-1 experiment (log(t+5.1), 4 data points, 128 virtual points). For each virtual point, estimate the fraction of the 50,000 RLRTO samples for which f'(s)=0 (i.e., the boundary is active). Then, using the same data and hyperparameters, compute the truncated-Gaussian posterior predictive mean and 95% CI at the test grid (via the exact conditional moments or a long NUTS run), and compare with the RLRTO predictive mean and CI. If the boundary fraction is non-negligible (say >5%) in a region where the true derivative is positive, or if the RLRTO CI is narrower than the truncated CI by more than 10%, the two methods are solving different inference problems, and the claimed 'comparable accuracy' is not a comparison of samplers for the same constrained posterior.","verdict_should_be":"UNCHANGED","load_bearing_attack":"The central claim requires that solving Eq. (36) repeatedly samples from a valid constrained posterior. I argue this premise is not met for the stated target. For a scalar linear-Gaussian model (x ~ N(0,1), y|x ~ N(x,1)), the RLRTO solution is x* = max(0, Z) with Z ~ N(mu_post, sigma^2_post). This censored distribution has a point mass at zero equal to Phi(-mu_post/sigma_post), whereas the truncated-Gaussian posterior used as baseline (Eq. 29) has no such atom. Thus RLRTO targets a different distribution. The paper concedes this in Sec. 3.4 ('cannot be written in a conventional Bayesian sense'), but still calls the output a 'constrained posterior' (Abstract) and compares predictive MSE/CI with methods targeting the truncated posterior. Moreover, the 'implicit prior' is data-dependent: the boundary atom's weight depends on the posterior covariance, so there may be no fixed prior for which Eq. (36) yields a posterior via Bayes' rule. Without evidence on the calibration of the RLRTO posterior or on how much it differs from the truncated target in the tested settings, the empirical 'comparable accuracy' could simply reflect that both methods enforce monotonicity, not that RLRTO is a valid sampler for the constrained posterior. The single-run metrics also prevent assessing the significance of the reported differences.","agreement_with_reader":"agree"},"referee_report":{"model":"deepseek-v4-flash","summary":"The paper develops a virtual-point Gaussian process framework for enforcing monotonicity constraints. It introduces a sampler based on regularized linear randomize-then-optimize (RLRTO), in which constrained derivative samples are obtained by repeatedly solving a non-negatively constrained randomized least-squares problem (Eq. 36). It also replaces the Gibbs samplers of two existing virtual-point methods (truncated prior and non-Gaussian likelihood) with NUTS. The methods are compared on six synthetic test functions and two differential-equation surrogate problems using MSE, 95% credible-interval width, integrated autocorrelation time, and effective samples per second. The authors report that RLRTO matches or improves predictive accuracy while being substantially faster, and they release code at a Zenodo repository.","tokens_in":19683,"tokens_out":5556,"duration_ms":88903,"significance":"The main contribution is the introduction of a fast, independent-sampling alternative to MCMC for constrained virtual-point GPs, with reproducible code and a broad experimental study. The RLRTO idea is interesting and the efficiency gains are large and credible. However, the statistical target of RLRTO is not characterized in the manuscript, and the sampling guarantee is imported from previous work; until the target distribution is specified and empirically validated, the comparison in Section 4 is not a comparison of samplers targeting the same posterior. The NUTS acceleration of the two baseline methods is a useful practical contribution, and the differential-equation demonstrations show the method's scope.","major_comments":[{"comment":"The statistical target of RLRTO is left implicit. The text states that the posterior 'can not be written in a conventional Bayesian sense', and the correctness of the sampler is cited to refs. [22–25], two of which are co-authored by an author of this manuscript. In a scalar linear-Gaussian example with a non-negativity constraint, the RLRTO solution has a censored distribution with a point mass at the boundary, whereas the truncated-Gaussian posterior in Eq. (29) has a continuous density. The paper therefore does not establish that RLRTO samples from the same posterior as the baselines; Section 4 and Tables 3–4 may be comparing different posterior models. Please provide a formal description of the distribution sampled by Eq. (36), including its dependence on the data and prior, or add a calibration/coverage experiment against a known target, or explicitly reframe the contribution as defining a new implicit-prior model rather than sampling an existing one.","section":"Sec. 3.4, Eq. (36)"},{"comment":"Equation (33) is not a valid joint Gaussian distribution because ReLU is a nonlinear transformation; the vector [f(t); ReLU(f'(s))] is not jointly Gaussian. The intended model appears to be f(t)|f'(s) ~ N(A ReLU(f'(s)), Σ*) with f'(s) ~ N(0, K11(s,s)), which yields the unnormalized density in Eq. (34). Please correct the displayed joint distribution and make the conditional likelihood explicit.","section":"Sec. 3.3, Eq. (33)"},{"comment":"All reported metrics appear to come from a single run without repeated trials or uncertainty estimates. Because the MSE and CI widths of the constrained methods are often close (e.g., the top-left panel of Fig. 5), the claims of 'comparable' or 'generally best' accuracy need support from multiple independent runs, different data realizations, or at least bootstrap confidence intervals. The efficiency comparisons are based on larger differences and are less sensitive, but the quality comparisons should be quantified.","section":"Sec. 4, Figs. 5–9 and Tables 3–4"}],"minor_comments":[{"comment":"References [22] and [44] are identical (Bardsley and Fox, 2012); please merge or renumber.","section":"References"},{"comment":"In Eq. (27), the covariance matrix is denoted K00(t,t) but K00 is never defined; this should presumably be K(t,t).","section":"Sec. 3.2, Eq. (27)"},{"comment":"The sentence 'based on the scaling behavior of the IAT seen in Fig. 6' appears to refer to Fig. 8, which displays IAT.","section":"Sec. 4.4"},{"comment":"The phrase 'Withe these samples' contains a typo.","section":"Sec. 4.1"},{"comment":"The warm-starting of the RLRTO solver, mentioned in the efficiency discussion in Section 4.4, is not described in the implementation details of Section 4.1; if warm starts are used, the resulting samples are not strictly independent and the IAT and ESS/s values should be interpreted accordingly.","section":"Sec. 4.4"}],"recommendation":"major_revision","confidential_remarks":"The main unresolved issue is the identity of the RLRTO target distribution; because two of the cited correctness guarantees are from the authors' own group, the revision should be held to a higher standard of proof or independent validation. The duplicate references [22]/[44] and the Eq. (33) error are easily fixed."},"author_rebuttal":null,"desk_editor":{"model":"deepseek-v4-flash","letter":"Worth a look. The genuinely new thing is the RLRTO/virtual-point combination: independent samples, no burn-in, and an order-of-magnitude improvement in ESS/s, backed by code on Zenodo. The NUTS substitution for Gibbs is less novel but sensible and well executed. I believe the efficiency results; the scaling plots are believable and the implementation uses standard trusted pieces (CUQIpy, Pyro, SciPy). The two differential-equation surrogates are a good stress test.\n\nThe soft spot is the one the stress test flags, and it does not go away on reading: RLRTO does not draw from the truncated-Gaussian posterior used by the truncated-prior baselines. In the scalar linear-Gaussian case the RLRTO solution is max(0, unconstrained posterior draw), so it puts point mass at zero with probability equal to the censoring probability; the truncated Gaussian has no such atom. The paper effectively concedes this in Section 3.4, where it says the posterior cannot be written in a conventional Bayesian sense, and Table 1 lists the density as implicit. But the abstract calls the output a constrained posterior and Section 4 compares MSE/CI width directly with methods that do target the truncated posterior. That is a comparison of three different posterior models, not three samplers of the same target. It does not invalidate the practical contribution—for surrogates built to enforce monotonicity, RLRTO clearly works—but the CI width comparisons are not apples-to-apples, and calibration of the RLRTO posterior is unexamined. A small test showing the RLRTO predictive distribution against the truncated-posterior distribution, or a calibration UQ plot, would settle it.\n\nTwo smaller issues. Equation (33) writes the joint distribution of f(t) and ReLU(f'(s)) as multivariate normal; that is false as stated since ReLU of a Gaussian is not Gaussian, and it is only the later density (34) that defines the model. The single-run metrics in Section 4 also make the small efficiency differences hard to assess; a few repeated data seeds would help.\n\nI don't worry about the self-references: the RLRTO guarantee is in prior peer-reviewed work, including one of the authors' own papers, but that is normal and the application here is independent. The paper deserves a serious referee, not a desk reject. I'd ask for a revision that clarifies the target distribution and adds a calibration check, then I'd accept it.","headline":"A genuinely useful and well-tested application of RLRTO to constrained GPs, whose practical speedup is real, but the paper should be clearer that RLRTO samples a different (censored) posterior than the truncated-Gaussian baselines.","tokens_in":20202,"tokens_out":3433,"would_cite":true,"duration_ms":35306,"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 claims that monotone Gaussian-process surrogates can be sampled by repeated constrained least-squares solves instead of MCMC, giving independent draws with autocorrelation near one.","keywords":["Gaussian processes","monotonicity constraints","virtual points","randomize-then-optimize","constrained posterior sampling","derivative-enhanced Gaussian processes","surrogate modeling","uncertainty quantification"],"falsifier":"Run a low-dimensional case, such as 1D-1 with four data points and a small virtual-point set, where exact samples from the truncated-Gaussian posterior are available by rejection or exact Monte Carlo; if a two-sample test on the derivative vectors shows the RLRTO empirical distribution does not approach the truncated-Gaussian target as the number of draws grows, then RLRTO samples a different posterior than intended.","tokens_in":19208,"feed_emoji":"📈","tokens_out":11624,"duration_ms":112702,"temperature":0.7,"pith_summary":"This paper claims that monotonicity constraints can be added to Gaussian-process surrogates much more cheaply than current MCMC practice by turning posterior sampling into repeated constrained least-squares solves. The proposed virtual-point method samples nonnegative derivative values at virtual points with the Regularized Linear Randomize-then-Optimize (RLRTO) scheme, which produces independent draws with no burn-in, no thinning, and no sequential correlation. On six synthetic functions and two differential-equation surrogates, it reports predictive accuracy comparable to or better than truncated-prior and non-Gaussian-likelihood baselines while sharply reducing sampling cost as the number of virtual points grows. The paper also shows that swapping Gibbs sampling for the No-U-Turn Sampler accelerates the two existing virtual-point methods. If the sampling premise holds, physically constrained GP surrogates become practical in settings where monotonicity is known in advance.","feed_headline":"Faster monotone GPs: each sample is one convex solve","feed_subtitle":"Independent draws with no burn-in beat Gibbs and NUTS on virtual-point Gaussian-process surrogates.","key_machinery":"The central object is RLRTO (Regularized Linear Randomize-then-Optimize), a sampling scheme that converts each posterior draw into one constrained quadratic program. The virtual-point construction places a finite set of constraint points $s$ in the input domain and enforces $f'(s)\\ge 0$ there; the derivative-enhanced GP supplies the linear forward map $A=K_{01}(t,s)K_{11}(s,s)^{-1}$ and the conditional formulas that convert a drawn derivative vector $f'(s)$ into predictive function draws. RLRTO perturbs the data vector with $\\hat b\\sim\\mathcal{N}(f(t),\\Sigma_*(t))$ and the prior mean with $\\hat c\\sim\\mathcal{N}(0,K_{11}(s,s))$, then solves the nonnegativity-constrained least-squares problem; the solution is the sample. The paper contrasts this with two MCMC alternatives, a truncated-Gaussian posterior sampled by Gibbs or NUTS and a ReLU-transformed non-Gaussian-likelihood posterior also sampled by Gibbs or NUTS.","core_discovery":"The paper's central claim is that the constrained posterior of a monotone Gaussian-process surrogate can be sampled, without MCMC, by repeatedly solving the constrained least-squares problem $\\min_{f'(s)\\ge 0}\\,\\frac{1}{2}\\|A f'(s)-\\hat b\\|_{\\Sigma_*^{-1}}^2+\\frac{1}{2}\\|f'(s)-\\hat c\\|_{K_{11}^{-1}}^2$, where $A=K_{01}(t,s)K_{11}(s,s)^{-1}$ links virtual-point derivative values to data values, $\\hat b$ is drawn from $\\mathcal{N}(f(t),\\Sigma_*(t))$, and $\\hat c$ is drawn from $\\mathcal{N}(0,K_{11}(s,s))$. Each solution is taken as one independent sample from the constrained posterior; the paper explicitly says this posterior is implicit and cannot be written in a conventional Bayesian closed form. The claim is that these optimization-based samples give predictive accuracy comparable to or better than the truncated-prior and non-Gaussian-likelihood baselines while making the sampling cost nearly independent of sample correlation, with integrated autocorrelation time near one and much higher effective samples per second. The same idea is reported to extend to other linear inequality constraints through projection mappings, and the paper also demonstrates that replacing Gibbs with the No-U-Turn Sampler accelerates the two baseline samplers.","pith_inferences":["The implicit-posterior caveat means that predictive accuracy parity with MCMC baselines does not by itself establish that RLRTO samples the truncated-Gaussian posterior; users who need faithful derivative-level credible intervals should first run a distributional check against a reference sampler.","If the optimization-based samples are accepted as the target, the practical ceiling is set by convex-solver speed rather than mixing time, which suggests monotone GP surrogates with thousands of virtual points become feasible in Bayesian optimization and PDE-inversion loops.","A natural testable extension is to enforce convexity in higher dimensions by replacing the scalar constraint $f''(s)\\ge 0$ with a positive-semidefinite Hessian constraint inside each randomized solve, a direction the paper leaves open."],"forward_implications":["Monotonicity-constrained GP surrogates can be built from independent samples, so no burn-in, thinning, or sequential correlation is needed, and sampling cost in effective samples per second improves sharply as the number of virtual points grows.","Flat regions are representable: unlike the truncated-prior method, RLRTO solutions can lie on the boundary $f'(s)=0$, matching the non-Gaussian-likelihood flexibility at lower cost.","The NUTS upgrade makes the two existing virtual-point methods practical for larger problems and, with enough virtual points, NUTS beats Gibbs on integrated autocorrelation time and effective samples per second.","Because RLRTO samples are independent, the sampling loop is embarrassingly parallelizable, so wall-clock time can be reduced further by running many randomized solves concurrently.","The constrained-optimization sampling extends to other linear inequality constraints (boundedness, convexity) and to higher-order derivative constraints, because only the projection onto the feasible set changes."],"supporting_citations":[{"why":"The original RLRTO formulation for nonnegativity-constrained posteriors, cited together with [23] for the constrained least-squares sampling problem.","marker":"[44]"},{"why":"Provides the RLRTO extension for nonnegativity-constrained linear inverse problems and the constrained least-squares formulation used as Eq. (36).","marker":"[23]"},{"why":"Supplies the projected-density interpretation used to justify treating constrained optimization solutions as posterior samples.","marker":"[24]"},{"why":"Develops sparse Bayesian inference with regularized Gaussian distributions, grounding the implicit-prior view of RLRTO solutions and the projection-mapping extension.","marker":"[25]"},{"why":"Introduces the Randomize-then-Optimize method for posterior sampling on which the paper's Algorithm 1 is built.","marker":"[42]"},{"why":"Establishes that RTO delivers independent samples in linear-Gaussian settings, supporting the efficiency claim for the GP application.","marker":"[43]"},{"why":"Introduces the virtual point method for monotonicity constraints that the framework refines.","marker":"[11]"},{"why":"Proposes the non-Gaussian likelihood model and its Gibbs sampler, one of the baselines the paper accelerates with NUTS.","marker":"[13]"},{"why":"Supplies the No-U-Turn Sampler used to replace Gibbs in the truncated-prior and non-Gaussian-likelihood baselines.","marker":"[26]"}],"fun_headline_variants":["Monotone GP sampling: one convex solve per draw","Skip MCMC: constrained GP via randomized solves","Independent GP samples without Gibbs or NUTS","Fast monotone surrogates: solve, don't sample","Constrained GP sampling as randomized least squares"],"cache_read_input_tokens":3200,"weakest_assumption_plain":"The load-bearing premise is that repeatedly solving the randomized constrained least-squares problem yields samples from the constrained posterior the user intends; the paper does not derive that posterior in closed form.","fun_headline_variants_meta":{"raw":{"variants":["Monotone GP sampling: one convex solve per draw","Skip MCMC: constrained GP via randomized solves","Independent GP samples without Gibbs or NUTS","Fast monotone surrogates: solve, don't sample","Constrained GP sampling as randomized least squares"]},"model":"deepseek-v4-flash","effort":"low","cost_usd":0.000214,"raw_usage":{"total_tokens":1474,"prompt_tokens":1044,"completion_tokens":430,"prompt_tokens_details":{"cached_tokens":384},"prompt_cache_hit_tokens":384,"prompt_cache_miss_tokens":660,"completion_tokens_details":{"reasoning_tokens":354}},"tokens_in":660,"tokens_out":430,"duration_ms":4521,"temperature":1.0,"reasoning_tokens":354,"cache_read_input_tokens":384,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-08-06T18:57:13.247088+00:00","model_set":{"reader":"deepseek-v4-flash"},"falsifier":"Run a low-dimensional case, such as 1D-1 with four data points and a small virtual-point set, where exact samples from the truncated-Gaussian posterior are available by rejection or exact Monte Carlo; if a two-sample test on the derivative vectors shows the RLRTO empirical distribution does not approach the truncated-Gaussian target as the number of draws grows, then RLRTO samples a different posterior than intended.","supporting_citations":[{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"The original RLRTO formulation for nonnegativity-constrained posteriors, cited together with [23] for the constrained least-squares sampling problem."},{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"Provides the RLRTO extension for nonnegativity-constrained linear inverse problems and the constrained least-squares formulation used as Eq. (36)."},{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"Supplies the projected-density interpretation used to justify treating constrained optimization solutions as posterior samples."},{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"Develops sparse Bayesian inference with regularized Gaussian distributions, grounding the implicit-prior view of RLRTO solutions and the projection-mapping extension."},{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"Introduces the Randomize-then-Optimize method for posterior sampling on which the paper's Algorithm 1 is built."},{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"Establishes that RTO delivers independent samples in linear-Gaussian settings, supporting the efficiency claim for the GP application."},{"cited_title":"Riihimaki, A","cited_arxiv_id":null,"evidence_quote":"Introduces the virtual point method for monotonicity constraints that the framework refines."},{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"Proposes the non-Gaussian likelihood model and its Gibbs sampler, one of the baselines the paper accelerates with NUTS."},{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"Supplies the No-U-Turn Sampler used to replace Gibbs in the truncated-prior and non-Gaussian-likelihood baselines."}],"review_version":1}