{"id":"bb789402-f9bc-4a1b-99cf-a24c82e319b0","arxiv_id":"2411.14570","paper_version":1,"verdict":"CONDITIONAL","confidence":"MODERATE","novelty_score":6.0,"correctness_risk":"medium","formal_verification":"none","parameter_count":3,"one_line_summary":"GradVI optimizes the variational empirical Bayes regression objective with quasi-Newton methods, achieving similar accuracy to CAVI but faster convergence on correlated predictors and much faster trend filtering.","lead":"This paper introduces GradVI, a gradient-based alternative to coordinate ascent for fitting variational empirical Bayes sparse regression models. It aims to make Bayesian sparse regression and trend filtering much faster by replacing one-at-a-time updates with quasi-Newton steps that rely mostly on matrix-vector products.","discovery_kind":"new_method","skeptic_critique":{"model":"deepseek-v4-flash","headline":"Paper's own experiments show the recommended Compound variant and the CAVI trend-filtering reference failing to converge within 2000 iterations, so the abstract's 'converges in fewer iterations when predictors are highly correlated' is overstated.","rationale":"I read the manuscript in good faith. The core idea is credible: the PLR reformulation of Kim et al. is used correctly, the Direct and Compound variants are principled ways to handle the intractable penalty, and the per-iteration complexity argument for structured X is coherent. The reader's weakest assumption—invertibility of S—turns out to be supportable: for normal-mixture and point-normal priors, the posterior mean in the normal means model is strictly increasing (monotone likelihood ratio), so T exists and the trisection inversion is defensible; the absence of proof is a presentation gap, not a correctness failure. The real load-bearing concern is empirical: the abstract's headline convergence claim conflicts with the paper's own reports of non-convergence for the recommended Compound variant in correlated settings, and the trend-filtering benchmark uses a non-converged CAVI run as reference. The Discussion's n=1e6 scalability remark is also unsupported. These issues do not invalidate the method, but they do mean the central claims should be qualified. The reader's CONDITIONAL verdict already captures this, so I leave the verdict unchanged.","tokens_in":16408,"tokens_out":12313,"duration_ms":124771,"concrete_test":"Using the public DSC repository (github.com/banskt/gradvi-experiments), extract per-replicate convergence flags, iteration counts, and final ELBOs for the center and right panels of Fig. 2. Compute the fraction of GradVI Compound and CAVI runs that hit the 2000-iteration cap, and the ELBO gap at cap. Then re-run the correlated block design with a 5000-iteration cap and a tolerance-based stopping rule (relative ELBO change < 1e-6) for both methods. If Compound's non-convergence rate is non-negligible while CAVI converges, or if its converged ELBO does not match the CAVI optimum, the abstract must be qualified to the convergent subset. Also run one n=1e6 trend-filtering benchmark to support the Discussion's scalability claim.","verdict_should_be":"UNCHANGED","load_bearing_attack":"The strongest claim in the Abstract and Sec. 1 is that GradVI converges in fewer iterations than CAVI when predictors are highly correlated and is much faster when X admits fast matrix-vector products. The paper's own results undercut this claim. In Sec. 4.1, for block-correlated predictors, 'GradVI Compound sometimes failed to converge after 2000 iterations', even though Compound is the variant recommended as default in Sec. 5. In Sec. 4.2, the trend-filtering comparison has no converged baseline: CAVI did not converge after 2000 iterations on the unscaled H matrix, and GradVI did not converge on the scaled H matrix; the ELBO 'reference' is a non-converged CAVI run. Thus the headline 'fewer iterations under high correlation' is only established for the subset of runs that terminate, and the speed comparison is confounded by comparing against a non-converged baseline. This matters because the practical case for GradVI rests on reliable convergence in correlated designs, not just on per-iteration complexity. The reader's inverse-T concern is less severe: for the ash and point-normal priors used here, S_g,v is strictly increasing by monotone likelihood ratio, so T exists and trisection is justifiable; the missing proof is a presentation gap. The more load-bearing issue is that the convergence evidence behind the central claim is incomplete.","agreement_with_reader":"partial"},"referee_report":{"model":"deepseek-v4-flash","summary":"The paper proposes GradVI, a quasi-Newton (L-BFGS-B) alternative to coordinate-ascent variational inference (CAVI) for fitting variational empirical Bayes sparse multiple regression with normal-mixture priors. It builds on the Kim et al. penalized-likelihood reformulation of the VEB objective, and develops two implementations: Direct (numerically inverting the posterior mean operator) and Compound (reparameterizing via the inverse operator). The paper reports simulation comparisons in high-dimensional regression and zeroth-order trend filtering, and claims that GradVI matches CAVI accuracy, converges in fewer iterations under correlated predictors, and is much faster when fast matrix-vector products are available. The central mathematical identities are the derivative formulas (22)-(24) and the compound objective (25).","tokens_in":16666,"tokens_out":5203,"duration_ms":47482,"significance":"If the claims held, GradVI would be a useful and flexible algorithmic contribution: it reduces the per-iteration cost to matrix-vector products, accommodates automatic differentiation and new priors, and the authors provide open-source code and reproducible DSC pipelines. The derivative identities follow from Tweedie's formula, and the PLR theorem is a published result from Kim et al., so the derivation is not circular. The practical significance is currently only partially established, however, because several empirical comparisons lack a converged baseline and the recommended variant sometimes fails to converge within the iteration budget.","major_comments":[{"comment":"The trend-filtering experiments do not establish the abstract's claim that GradVI 'converges in fewer iterations' or 'is much faster than CAVI' in settings with fast matrix-vector products. The text states that CAVI did not converge within 2000 iterations on the unscaled H matrix, that GradVI Compound did not converge within 2000 iterations on the scaled H matrix, and that the ELBO/RMSE 'reference' is a non-converged CAVI run. The scaled-H comparison also changes the design matrix and hence the induced prior on the regression coefficients, so ELBO values across scaled and unscaled runs are not directly comparable. Please report convergence rates for every method/design, compare against a converged baseline (for example, a longer CAVI run or a smaller ELBO tolerance), or explicitly qualify the speed and iteration claims as applying only to runs that terminate.","section":"Sec. 4.2, Fig. 2 (right column); Fig. 5"},{"comment":"The recommendation of Compound as the default is undercut by its failure to converge in some correlated-predictor simulations. The text reports 'GradVI Compound sometimes failed to converge after 2000 iterations, whereas GradVI Direct converged in all simulations.' Because the abstract's 'converges in fewer iterations when the predictors are highly correlated' claim is based on these simulations, the paper must report how often each variant fails, the distribution of iteration counts conditional on convergence, and whether the iteration summaries in Fig. 2 include non-converged runs (for example, censored at 2000). Without this, the claim holds only for a selected subset of runs.","section":"Sec. 4.1 (Block-correlated variables); Sec. 5"},{"comment":"The existence and properties of the inverse operator T = S^{-1} are assumed rather than proved for the general case. Symmetry and monotonicity are asserted only for constant v_j, yet the Direct method applies trisection for coordinate-dependent v_j and the Compound reparameterization z = T(theta) requires T to be a bijection over the relevant range. Please provide a proof (or precise conditions) that S_{g,v_j} is strictly increasing and onto for the ash and point-normal priors used here, and give bracketing/error guarantees for the trisection inversion. If such a guarantee is not available, state the resulting limitation explicitly.","section":"Sec. 3.1 and 3.2, Eqs. (19)-(25)"}],"minor_comments":[{"comment":"There are several typos, including 'optimzation' in Sec. 3.1 and 'tractible' in Secs. 3.2 and 5; please correct them.","section":"Sec. 3.1 and Sec. 3.2"},{"comment":"The notation paragraph says vectors are denoted by bold lowercase letters with 'e.g., b', but then writes 'e.g., y or Y'; please make the vector/scalar notation for y consistent.","section":"Sec. 2.1"},{"comment":"The text refers to panels (a), (b), and (c) of Fig. 4, but the figure as described does not show clear panel labels; please add labels to the figure panels.","section":"Sec. 4.2 and Fig. 4"},{"comment":"The Discussion claims GradVI 'enables trend filtering for datasets with n = 10^6 on a single core with 100 GB of memory,' but no experiment in the manuscript supports this; please indicate whether this is a reported run or an extrapolation.","section":"Sec. 5"}],"recommendation":"major_revision","confidential_remarks":"The manuscript is honest about its non-convergence issues in the text, but the abstract and introduction overstate the convergence and speed advantages relative to the evidence. The missing convergence rates and converged baseline are fixable within the manuscript's scope, so I see this as a major-revision rather than a reject. The relationship to Kim et al. is properly cited; the methodological novelty is moderate but the software and reproducible workflow are strengths."},"author_rebuttal":null,"desk_editor":{"model":"deepseek-v4-flash","letter":"First thing you should know: this is a genuinely useful algorithmic paper, not a new statistical model. The idea is to replace coordinate ascent with L-BFGS-B on the penalized-regression objective from Kim et al. for variational empirical Bayes sparse regression. The two tricks for dealing with the non-closed-form penalty—numerical inversion of the posterior-mean operator, and a compound-penalty reparametrization that avoids inversion—are sound. The derivative identities (22)-(24) check out via Tweedie's formula. The code and reproducibility scaffolding are real, and I believe the method will be used.\n\nThe main new thing is the compound reparametrization plus the use of fast matrix-vector products. The paper shows convincingly that in settings where X admits a fast matvec, GradVI is much faster per iteration than CAVI, and that it often reaches a better ELBO. The trend-filtering experiments are a good demonstration, and the comparison with genlasso is a nice addition. Credit where due: the paper is honest about some of its own failures, like the Compound variant not always converging in the correlated-predictor setting.\n\nWhere I'd push back is the abstract's headline claim that GradVI 'converges in fewer iterations when the predictors are highly correlated.' Their own experiments undercut that. In the block-correlated setting, the Compound variant—the one they recommend as default—sometimes failed to converge within 2000 iterations, and the Direct variant was the one that always converged. More importantly, in the trend-filtering compare, the CAVI reference did not converge on the unscaled H matrix, so the speed and ELBO comparisons are against a non-converged baseline. That doesn't kill the paper—the method still has clear value for fast-matvec problems and for less initialization sensitivity—but it means the abstract overstates the convergence advantage. The claim should be conditioned on the settings where it actually holds.\n\nThe inverse-T concern the stress-test flagged is less severe. For the ash and point-normal priors used, the posterior-mean operator is monotone by monotone likelihood ratio, so trisection is justifiable. The missing proof is a presentation gap, not a flaw in the algorithm.\n\nBottom line: this is a methods paper with solid experiments but overstated convergence claims. It deserves a serious referee; with revision to temper the abstract and add a convergence table that distinguishes converged vs not, it'd be a useful addition to the optimization-for-VI literature. I'd cite it for the compound-penalty trick. Not the strongest paper of the year, but honest, reproducible, and practical.","headline":"A practical, honest methods paper on gradient-based VEB regression; the compound-penalty trick is real, but the convergence claims are overreaching.","tokens_in":17215,"tokens_out":2742,"would_cite":true,"duration_ms":22367,"reading_group":"yes","serious_thinker":"yes","would_accept_peer_review":true},"rs_alignment":null,"lean_confirmation":null,"pith_extraction":{"msc":["62F15","62J07","62J05","65K05"],"pacs":[],"model":"deepseek-v4-flash","headline":"Replacing coordinate ascent with quasi-Newton steps keeps predictive accuracy, cuts iterations on correlated predictors, and reduces each optimization step to matrix-vector products.","keywords":["variational inference","empirical Bayes","multiple regression","penalized linear regression","quasi-Newton optimization","trend filtering","normal means model","posterior mean operator"],"falsifier":"Take a spike-and-slab prior and compute $S_{g,v}(z)$ on a fine grid of $z$ values for several variances $v$; if any observed $\\theta$ has more than one preimage under $S$, then the inverse $T$ used by both variants is not single-valued, and Direct's numerical inversion or Compound's reparameterization would fail on a regression with that prior.","tokens_in":1935,"feed_emoji":"📈","tokens_out":3634,"duration_ms":113925,"temperature":0.7,"pith_summary":"This paper proposes gradient-based variational inference (GradVI) for sparse multiple regression, replacing the coordinate-ascent updates that most variational empirical Bayes (VEB) methods use. The central claim is that the VEB objective can be minimized with a quasi-Newton method, achieving predictive accuracy similar to coordinate ascent while using fewer iterations when predictors are correlated. Because the per-iteration work reduces to matrix-vector products, GradVI becomes much faster than coordinate ascent whenever the design matrix admits fast multiplication, such as in trend filtering. The paper develops two variants to handle the fact that the penalty in the penalized-regression view is not available in closed form: one inverts the normal-means posterior mean operator numerically, and the other reparameterizes to make the penalty analytic.","feed_headline":"Gradient updates speed up Bayesian regression on correlated data","feed_subtitle":"Replacing coordinate ascent with quasi-Newton steps retains accuracy and unlocks fast trend filtering.","key_machinery":"The central object is the normal-means posterior-mean operator $S_{g,v}(z)$, which maps a noisy observation $z$ to the posterior mean of the signal under prior $g$ and variance $v^2$; its inverse is denoted $T$. The penalized-regression theorem converts VEB into minimizing $\\frac{1}{2\\sigma^2}\\|y - X\\theta\\|^2$ plus a sum of penalties $\\rho(\\theta_j, g, v_j^2)$, where the penalty satisfies $\\rho'(\\theta_j) = (T(\\theta_j) - \\theta_j)/v_j^2$. This identity means that the intractable penalty and its gradients can be evaluated from $T$ and the normal-means log marginal likelihood $\\ell_{NM}$. The paper's two algorithms are two ways to compute $T$: Direct uses trisection-style numerical inversion, with a faster switch-and-spline inversion when $v_j$ is constant, and Compound reparameterizes through $z = T(\\theta)$ to obtain an analytic penalty. A limited-memory quasi-Newton optimizer then minimizes either objective, so the dominant per-step work is the matrix-vector products $Xr$ and $X^T r$.","core_discovery":"The paper builds on the result that variational empirical Bayes regression with a factorized posterior is exactly a penalized regression in the posterior mean $\\theta$: minimizing $h(\\theta, g, \\sigma^2) = \\frac{1}{2\\sigma^2}\\|y - X\\theta\\|^2 + \\sum_j \\rho(\\theta_j, g, v_j^2) + \\text{constants}$, where the penalty $\\rho$ is defined through the normal-means posterior mean $S_{g,v}(z)$ and its inverse $T = S^{-1}$. The discovery is that this penalized regression objective can be optimized by quasi-Newton methods rather than coordinate updates, provided one can evaluate $T$. The Direct variant evaluates $T$ by numerical inversion at each evaluation; the Compound variant changes variables to $z = T(\\theta)$ so that $\\theta = S(z)$, making the compound penalty analytically tractable and avoiding inversion during optimization. The paper shows numerically that both variants match coordinate-ascent predictive accuracy, that the Compound variant often reaches a better ELBO in fewer iterations on block-correlated designs, and that in trend filtering the gradient methods complete far more iterations within the same time because the matrix-vector products $Xr$ and $X^T r$ exploit the structure of the design matrix.","pith_inferences":["The same penalized-regression reformulation should extend beyond normal means: any variational family whose mean map has a tractable inverse would let GradVI-style optimization handle non-Gaussian likelihoods, such as logistic regression, a natural extension the paper does not test.","Because the cost is dominated by matrix-vector products, GradVI could be combined with randomized linear algebra, sketching, or fast transforms to attack designs that are too large to store explicitly.","The Compound variant's reliance on a bijective map $S$ suggests a practical diagnostic: before using it with a new prior, one should check monotonicity of $S$ over the range of $v_j$ values that arise in optimization, since the paper only gives a monotonicity argument for constant $v_j$.","The comparison with cross-validated $\\ell_1$ trend filtering hints at a general phenomenon: empirical-Bayes penalties that estimate component sizes from data may recover signal magnitudes that a single $\\ell_1$ penalty systematically shrinks."],"forward_implications":["Any design matrix with a fast multiplication routine, such as trend filtering, Fourier or wavelet bases, or sparse graphs, becomes a setting where GradVI is much faster than coordinate ascent, because both objective and gradient evaluations reduce to $Xr$ and $X^T r$ products.","For correlated predictors, GradVI with a null initialization can match or beat coordinate ascent initialized from a lasso fit, which removes the costly cross-validation step that large regression problems currently require.","The modular treatment of priors means that adding a new prior family only requires solving the corresponding normal-means problem, rather than deriving and debugging full coordinate updates for that prior.","In Bayesian trend filtering, the paper reports better estimates of changepoint magnitudes than cross-validated $\\ell_1$ trend filtering, especially at high noise, and reports scalability to $n = 10^6$ observations on a single machine.","The gradient-based structure can exploit parallel and GPU-accelerated matrix multiplication in ways that sequential coordinate ascent cannot."],"supporting_citations":[{"why":"Supplies the penalized-regression theorem that rewrites the VEB objective as a penalized regression with penalty defined through the normal-means posterior mean; it is the foundation of both GradVI variants and the source of the coordinate-ascent baseline.","marker":"[29]"},{"why":"Defines the adaptive shrinkage (ash) normal-mixture prior used in all numerical experiments and as the default prior module in the software.","marker":"[28]"},{"why":"Introduces $\\ell_1$ trend filtering, the nonparametric setting whose structure provides the fast-matrix-vector test bed for GradVI.","marker":"[4]"},{"why":"Proves the equivalence of trend filtering to a lasso problem with the matrix $H^{(k+1)}$, which is what allows GradVI to exploit fast matrix-vector products.","marker":"[5]"},{"why":"Supplies the mean-field variational approximation and coordinate-ascent variational inference framework that GradVI is compared against.","marker":"[12]"},{"why":"Tweedie's formula gives the identity $S(z) = z + s^2 \\ell'_{NM}(z)$, used to derive the penalty derivatives.","marker":"[41]"},{"why":"Supplies the limited-memory BFGS optimization routine used to minimize the penalized regression objectives.","marker":"[44]"}],"fun_headline_variants":["Quasi-Newton steps speed up variational Bayes regression","GradVI: quasi-Newton beats coordinate ascent for Bayes regression","Gradient-based variational inference accelerates Bayesian regression","Quasi-Newton optimization gives faster Bayesian regression on correlated data","GradVI: gradient-based Bayesian regression that scales with structure"],"cache_read_input_tokens":19328,"weakest_assumption_plain":"Everything rests on the assumption that the normal-means posterior-mean map $S$ can be inverted reliably on every value the optimizer visits: Direct numerically solves $S(z) = \\theta$, and Compound assumes $\\theta \\mapsto S^{-1}(\\theta)$ is a well-defined bijection, an assumption that is justified only for constant variances and otherwise relies on bracketing and monotonicity that the paper does not verify.","fun_headline_variants_meta":{"raw":{"variants":["Quasi-Newton steps speed up variational Bayes regression","GradVI: quasi-Newton beats coordinate ascent for Bayes regression","Gradient-based variational inference accelerates Bayesian regression","Quasi-Newton optimization gives faster Bayesian regression on correlated data","GradVI: gradient-based Bayesian regression that scales with structure"]},"model":"deepseek-v4-flash","effort":"low","cost_usd":0.000799,"raw_usage":{"total_tokens":3572,"prompt_tokens":1058,"completion_tokens":2514,"prompt_tokens_details":{"cached_tokens":384},"prompt_cache_hit_tokens":384,"prompt_cache_miss_tokens":674,"completion_tokens_details":{"reasoning_tokens":2433}},"tokens_in":674,"tokens_out":2514,"duration_ms":15526,"temperature":1.0,"reasoning_tokens":2433,"cache_read_input_tokens":384,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-08-12T15:08:45.162758+00:00","model_set":{"reader":"deepseek-v4-flash"},"falsifier":"Take a spike-and-slab prior and compute $S_{g,v}(z)$ on a fine grid of $z$ values for several variances $v$; if any observed $\\theta$ has more than one preimage under $S$, then the inverse $T$ used by both variants is not single-valued, and Direct's numerical inversion or Compound's reparameterization would fail on a regression with that prior.","supporting_citations":[{"cited_title":"A flexible empirical Bayes approach to multiple linear regression and connections with penalized regression","cited_arxiv_id":null,"evidence_quote":"Supplies the penalized-regression theorem that rewrites the VEB objective as a penalized regression with penalty defined through the normal-means posterior mean; it is the foundation of both GradVI variants and the source of the coordinate-ascent baseline."},{"cited_title":"False discovery rates: a new deal","cited_arxiv_id":null,"evidence_quote":"Defines the adaptive shrinkage (ash) normal-mixture prior used in all numerical experiments and as the default prior module in the software."},{"cited_title":"ℓ1 Trend Filtering","cited_arxiv_id":null,"evidence_quote":"Introduces $\\ell_1$ trend filtering, the nonparametric setting whose structure provides the fast-matrix-vector test bed for GradVI."},{"cited_title":"Adaptive piecewise polynomial estimation via trend filtering","cited_arxiv_id":null,"evidence_quote":"Proves the equivalence of trend filtering to a lasso problem with the matrix $H^{(k+1)}$, which is what allows GradVI to exploit fast matrix-vector products."},{"cited_title":"Scalable Variational Inference for Bayesian Variable Selection in Regression, and Its Accuracy in Genetic Association Studies","cited_arxiv_id":null,"evidence_quote":"Supplies the mean-field variational approximation and coordinate-ascent variational inference framework that GradVI is compared against."},{"cited_title":"Proceedings of the Third Berkeley Symposium on Mathematical Statistics and Probability, vol. 1: Contributions to the Theory of Statistics","cited_arxiv_id":null,"evidence_quote":"Tweedie's formula gives the identity $S(z) = z + s^2 \\ell'_{NM}(z)$, used to derive the penalty derivatives."}],"review_version":1}