Pith. sign in

REVIEW 3 major objections 6 minor 1 cited by

hdbayes: An R Package for Bayesian Analysis of Generalized Linear Models Using Historical Data

T0 review · 3 major / 6 minor · reviewed 2026-08-06 · deepseek-v4-flash

Pith's one-line read The paper introduces hdbayes, an R package that implements seven historical-data borrowing priors for generalized linear models behind one consistent glm-style syntax, with Stan-based sampling and marginal-likelihood tools.

desk verdict A genuinely useful R package for historical-data priors, but the multi-historical NPP formula is structurally not a normalized power prior and should be fixed before publication. read the letter →

arxiv 2506.20060 v1 pith:WN4IMZZZ submitted 2025-06-24 stat.ME stat.CO

classification stat.MEstat.CO MSC 62F1562J12
keywords historicaldatapowerpriornormalizedcommensuratelatentexchangeabilitymeta-analyticpredictivegeneralizedlinearmodelsRpackage
verification ladder T0 review T1 audit T2 compute T3 formal

The pith

A machine-rendered reading of the paper's core claim, the machinery that carries it, and where it could break.

The reading

This paper claims that the main barrier to using historical-data priors in Bayesian regression is software fragmentation, and that a single R package can remove it. hdbayes implements seven such priors—the power prior, normalized power prior, normalized asymptotic power prior, Bayesian hierarchical model, robust meta-analytic predictive prior, commensurate prior, and latent exchangeability prior—for generalized linear models, all with essentially the same glm-style syntax and all sampling done by Stan's Hamiltonian Monte Carlo. The authors demonstrate the package on two clinical trial datasets, including a logistic regression comparison across all priors and a proportional-hazards analysis via the Poisson likelihood representation. If the package works as claimed, applied statisticians can try many borrowing strategies in one workflow and compare them fairly.

What carries the argument

The load-bearing object is the wrapper function glm.prior(formula, family, data.list, prior.args, ...), which mimics R's glm() and turns any of the seven priors into a Stan program. The computational workhorses are Hamiltonian Monte Carlo/NUTS sampling from Stan and bridge-sampling estimation of normalizing constants, with the grid-plus-LOESS-plus-linear-interpolation approximation used for the normalized power prior. The normalized asymptotic power prior uses the asymptotic normal approximation to the power prior under a log transformation of the dispersion parameter, giving a closed-form prior, while the latent exchangeability prior uses a finite mixture model to discount historical observations individually.

What would settle it

Fit the normalized power prior to a small normal linear model where the normalizing constant Z(a0) is known in closed form, using glm.npp with the package's grid estimates, and compare the resulting posterior for the discounting parameter and coefficients with the exact posterior; any systematic discrepancy that grows as the sample shrinks would show the approximation chain is unreliable.

Watch

Extended reading notes

Core claim

The central claim is that a cohesive implementation is possible: every one of these priors can be expressed as a small wrapper around a common formula, family, and data-list interface, with sensible defaults, and inference can be carried out in Stan without hand-tuning. The package's key methodological move is to handle the unnormalized priors by estimating normalizing constants with bridge sampling: for the normalized power prior, it samples the power-prior posterior on a grid of discounting values, smooths the estimated log-normalizing constants with LOESS, and feeds the smoothed curve into Stan as a lookup grid with linear interpolation; for the robust MAP prior, it estimates the marginal likelihoods of the informative and vague components and reweights their posterior samples. The paper also introduces a normalized asymptotic power prior that avoids normalizing-constant estimation entirely by using the asymptotic multivariate normal form of the power prior on log-dispersion.

Load-bearing premise

The load-bearing premise is that the package's numerical estimates of the constants that normalize the priors—produced by sampling on a grid of discount values, smoothing, and interpolating—remain accurate enough in real problems, since the paper demonstrates this only on a single dataset.

Editorial extensions

If this is right

  • Users can fit all seven priors with one consistent interface, so method comparisons no longer require translating between different packages and sampling schemes.
  • Marginal likelihoods for unnormalized priors become computable, enabling Bayes-factor model and link selection that was previously awkward.
  • Time-to-event analyses can be conducted through the Poisson representation of piecewise-constant-hazard proportional hazards models, extending the package beyond GLM responses.
  • Stan's NUTS sampling removes the need to tune Metropolis or slice samplers, reducing a common obstacle to using these priors in high dimensions.
  • The package can stack multiple historical data sets for power priors and the latent exchangeability prior, allowing a single analysis to pool several sources.

Reading between the lines

Editorial extensions of the paper, not claims the author makes directly.

  • The paper does not test the accuracy of the bridge-sampling-plus-smoothing normalizing constant estimates beyond one dataset, so users should check stability across grid choices and MCMC seeds before relying on normalized power prior or robust MAP results.
  • The same wrapper architecture could be extended natively to longitudinal or time-to-event models, and to propensity-score integrated priors, as the discussion notes.
  • The consistent interface could serve as a benchmark for future software: any new historical-data prior becomes adoptable if it fits the same formula, family, and data-list contract.
  • The empirical comparison on the AIDS data suggests that borrowing strength varies by method, so reporting several priors rather than one may be the honest default in clinical analyses.
Share X Bluesky LinkedIn Reddit HN

Signed reviews

No signed human review yet.

Editorial analysis

A structured set of objections, weighed in public.

Desk editor's note, referee report, and a circularity audit.

Referee Report

3 major / 6 minor

Summary. The paper introduces hdbayes, an R package that implements several historical-data priors for Bayesian generalized linear models: the power prior, normalized power prior, normalized asymptotic power prior, Bayesian hierarchical model, robust meta-analytic predictive prior, commensurate prior, and latent exchangeability prior. The package is built on Stan and cmdstanr with a consistent formula/family interface, and includes functionality for marginal-likelihood-based model selection. The manuscript reviews each prior, illustrates the package on an HIV/AIDS clinical trial dataset, and shows how time-to-event data can be handled through Poisson GLM representation. The paper also proposes numerical approaches for normalizing constants and mixture-weight updating, and provides code examples and output summaries.

Significance. If the implementations are correct, hdbayes is a valuable contribution: it is a CRAN-released package that unifies several historical-data priors for GLMs under a single, user-friendly, HMC-based interface, and it provides the first publicly available implementation of the latent exchangeability prior. The consistent syntax and inclusion of marginal likelihood tools are genuine strengths. However, the manuscript contains a load-bearing mathematical error in the definition of the multi-historical normalized power prior (Eq. 6), and it overstates the package's GLM coverage in Section 5 and Table 1. The numerical approximation chain for normalizing constants is also not validated beyond a single example, which is insufficient for a software paper claiming to implement these priors. These issues affect the reliability of the package for multi-historical settings and require correction before the paper can be considered for publication.

major comments (3)
  1. [Section 2.3, Eq. (6) and Algorithm 1] For H>1, the product of per-dataset normalized power priors in Eq. (6) is not a properly normalized density over (β, φ) conditional on a0. The integral of the product is C(a0) = ∫ π0(β,φ) ∏_{h=1}^H L_h(β,φ|D0h)^{a0h} dβ dφ / ∏_{h=1}^H Z_h(a0h), which generally depends on a0 and is not equal to 1. Consequently, the posterior over a0 under Eq. (6) is weighted by the omitted factor Z_joint(a0)/∏_h Z_h(a0h), where Z_joint(a0) is the joint normalizing constant. As a concrete counterexample, for two identical historical datasets the likelihood depends only on the total discount s = a01 + a02, but the product-normalized form gives different prior weight to splits such as (0.5,0.5) versus (0.1,0.9). Thus glm.npp for H>1 does not implement the normalized power prior as stated. Algorithm 1 must either be revised to use the joint normalizing constant or the package/paper should restrict the NPP to a single historical dataset and clearly disclose this limitation.
  2. [Section 2.2 vs Section 5 and Table 1] The paper claims in Section 5 that 'the full suite of GLMs and link functions available in the stats package is available in hdbayes' and Table 1 marks 'All models in stats::glm' for hdbayes. However, Section 2.2 states that binomial models with more than one trial are not supported by glm.pp. Since binomial models with a trial count greater than one are a standard part of stats::glm, these statements are contradictory and the 'full suite' claim is false as written. The paper should either document the actual supported family/link combinations precisely, or add support for binomial trials; at minimum, the overclaim in Section 5 should be removed.
  3. [Sections 2.3 and 4.1] The accuracy of the normalizing constant estimation chain for the NPP—bridge sampling on a grid of power-prior samples, LOESS smoothing, and linear interpolation inside Stan—is not validated beyond the single HIV example. The paper notes in Section 2.3 that for the normal linear model with a conjugate normal-gamma initial prior the normalizing constant is known, and lm.npp avoids estimation, but no comparison between glm.npp's bridge-sampling estimate and the analytic value is reported. No simulation study or independent check is provided to show that the estimated Z_h(a0h), and hence the posterior of a0 and β, is reliable. Since the NPP posterior depends critically on Z_h, the paper should include at least one validation study, such as recovery of known normalizing constants in the normal linear model or a comparison with an alternative implementation, to support the correctness of the approximation.
minor comments (6)
  1. [Section 2.2, paragraph on binomial models] The sentence 'The implementation allows for binomial models with the number of trials exceeding 1, which is a limitation of the implementation in the hdbayes implementation glm.pp' is garbled. It appears that 'allows' should be 'does not allow' or the sentence should be restructured to state clearly that glm.pp does not support binomial models with more than one trial while other implementations do.
  2. [Section 2.4, Eq. (7)] The expression for the asymptotic covariance matrix is typeset in a confusing way: 'a−1 0h h I(ˆβ0h, ˆϕ0h|D0h) i−1' should be written as a0h^{-1} [I(ˆβ0h, ˆϕ0h|D0h)]^{-1} to indicate the inverse of the Fisher information multiplied by the reciprocal of the discounting parameter.
  3. [Algorithm 1, line 4] The text says the bridge sampling estimate uses 'samples from the prior,' but steps 2–3 sample from the power prior (4) with the given discounting parameter. Please rephrase to 'samples from the power prior posterior' or similar for clarity.
  4. [Section 2.6, Eq. (14) and Algorithm 2] In Eq. (14), the denominator uses Z_V(D0), but the vague prior normalizing constant for the current data should be Z_V(D). Similarly, Algorithm 2 line 2 writes p_V(θ|D0) where the current-data posterior is intended. Also, the vague prior π_V is not explicitly defined in Section 2.6; please state that it is the same initial prior as in (5) or provide the formula.
  5. [Section 4.1, code description] The sentence 'In hdbayes, we assume that the first specified data set in the list is the current data sets' contains a typo: 'data sets' should be 'data set'.
  6. [Section 3, Eq. (18) and surrounding text] The description of the two-step marginal likelihood computation has awkward phrasing: 'after taking MCMC samples of the prior' and 'estimating the normalizing constant of the posterior' should be clarified to indicate that the MCMC samples are from the prior and posterior densities respectively, with bridge sampling applied to each.

Circularity Check

0 steps flagged · score 1.0 of 10

No significant circularity; the paper is a software-implementation report whose methods are defined from first principles or cited to prior external work, and whose self-citations are attributional rather than load-bearing.

full rationale

The paper's central content is an R package that implements several well-established historical-data priors (power prior, normalized power prior, Bayesian hierarchical model, robust MAP prior, commensurate prior, LEAP). The derivations shown are standard definitions or direct extensions of existing ones. The power prior in Eq. (3)-(4) is a textbook definition. The normalized power prior in Eq. (6) is defined explicitly, and Algorithm 1 implements that definition using per-dataset bridge-sampling estimates of the normalizing constants; the final posterior is computed by substituting the estimated functions into the given expression. This is internal to the method (estimating an intractable normalization), not a reuse of the target result. The RMAP procedure in Algorithm 2 computes mixture weights from bridge-sampling estimates of marginal likelihoods and then samples from the two component posteriors; this is a Monte Carlo composition, not a fitted parameter renamed as a prediction. Self-citations (LEAP to Alt et al. 2024a; NPP two-step approach to Carvalho and Ibrahim 2021; hdbayes package to Alt et al. 2024b) are credits for prior method development and do not carry the derivation. One mathematical concern is noted: for multiple historical datasets, Eq. (6) normalizes each per-dataset power prior separately rather than using the joint normalizing constant, so the resulting product is not a correctly normalized multi-historical NPP. However, this is a correctness issue about the definition and its implementation, not circularity: the implementation is faithful to the paper's own stated definition, and no fitted parameter or self-citation is masquerading as an independent result. Therefore the paper's claims do not reduce to their inputs by construction.

Assumptions & free parameters 3 free parameters · 3 assumptions · 0 invented entities

The paper introduces no new statistical entities: the priors, parameters, and hyperparameters all already exist in the prior literature. The package is a unification, so the ledger contains the user-chosen discounting and mixture weights and the approximation axioms behind the NPP and NAPP.

free parameters (3)
  • Discounting parameter a0 in the power prior = 0.11 in the HIV example
    In the HIV example, a0 is set to (1/2)(n/n0) as a heuristic, and the authors recommend trying several values. It is not fitted to the posterior target, but it is a user-choice input in the methods.
  • Beta prior parameters for a0 in NPP/NAPP = alpha=0.77, beta=6.21
    These are chosen to match a mean of (1/2)(n/n0) and a coefficient of variation of 1; an ad hoc elicitation, but not fitted to the result.
  • RMAP mixture weight gamma = 0.1 default
    The default is set in the package and the authors encourage testing other values; it is an elicited hyperparameter.
assumptions (3)
  • domain assumption The GLM likelihood as the data-generating model.
    Section 2, Eq. 1. All priors are defined with respect to this likelihood.
  • domain assumption For the NPP, the normalizing constant can be estimated by the two-step bridge sampling with LOESS smoothing and linear interpolation.
    Section 2.3, Algorithm 1. The accuracy of the smoothed function is demonstrated on one data example and is a load-bearing approximation.
  • domain assumption For the NAPP, the asymptotic normal approximation of the power prior is adequate.
    Section 2.4, Eq. 7. The paper explicitly states that the approximation may be overly informative when the normal is a poor fit.

how reviews work

0 comments
Cite this review

Pith. "Pith review of hdbayes: An R Package for Bayesian Analysis of Generalized Linear Models Using Historical Data." pith.science (2026). https://pith.science/paper/WN4IMZZZ

@misc{pith2026250620060,
  author       = {Pith},
  title        = {Pith review of: hdbayes: An R Package for Bayesian Analysis of Generalized Linear Models Using Historical Data},
  year         = {2026},
  howpublished = {\url{https://pith.science/paper/WN4IMZZZ}},
  note         = {Machine review of arXiv:2506.20060}
}
read the original abstract

There has been increased interest in the use of historical data to formulate informative priors in regression models. While many such priors for incorporating historical data have been proposed, adoption is limited due to access to software. Where software does exist, the implementations between different methods could be vastly different, making comparisons between methods difficult. In this paper, we introduce the R package hdbayes, an implementation of the power prior, normalized power prior, Bayesian hierarchical model, robust meta-analytic prior, commensurate prior, and latent exchangeability prior for generalized linear models. The bulk of the package is written in the Stan programming language, with user-friendly R wrapper functions to call samplers.

Figures

Figures reproduced from arXiv: 2506.20060 by the authors.

Figure 1
Figure 1. Beta prior distribution for a0 and estimated log￾normalizing constants as a function of a01 for the AIDS progression logistic regression model. In panel A), the vertical dotted line marks the value 1/2(n/n0). In panel B) we show the log-normalizing constant esti￾mated using Algorithm 1. Note that log Zh(0) = 0 for all h. As discussed in Section 2.3, for the NPP it is also necessary to esti￾mate the normalizing const… view at source ↗
Figure 2
Figure 2. Coefficient estimates for the AIDS progression logistic regression model. For each coefficient (vertical tiles) we show the posterior mean (dot) and 95% equal-tailed credible interval (solid vertical lines) ac￾cording to each of the models used here. Horizontal lines (dotted: ACTG036, dashed: ACTG019) mark the MLE for each coefficient. 4.2 Time-to-event example Because hdbayes focuses on the class of GLMs, it cannot… view at source ↗

Discussion (0). Continue with ORCID to comment.

Forward citations

Cited by 1 Pith paper

Reviewed papers in the Pith corpus that reference this work. Sorted by Pith novelty score. Full citation record

  1. Bayesian weighted discrete-time dynamic models for association football prediction

    stat.ME 2025-08 unverdicted novelty 6.0 of 10

    A Bayesian football goal model uses period-specific spike-and-slab priors to adaptively weight how much team attack and defense strengths change over time, improving prediction relative to standard dynamic models.

Reference graph

Works this paper leans on

15 extracted references · 7 canonical work pages · cited by 1 Pith paper

  1. [1]

    E. M. Alt, X. Chang, X. Jiang, Q. Liu, M. Mo, H. A. Xia, and J. G. Ibrahim. LEAP: The latent exchangeability prior for borrowing information from historical data.Biometrics, 80(3):ujae083, Sept. 2024a. ISSN 0006-341X. doi: 10.1093/biomtc/ujae083. E. M. Alt, X. Chen, L. M. Carvalho, and J. G. Ibrahim.hdbayes: Bayesian Analysis of Generalized Linear Models ...

  2. [1982]

    URLhttps://doi.org/ 10.1214/aos/1176345693. J. Gabry, R. Češnovar, A. Johnson, and S. Bronder.cmdstanr: R Interface to ‘CmdStan’,

  3. [1986]

    1986.10478240

    URLhttps://doi.org/10.1080/01621459. 1986.10478240. P. A. Volberding, S. W. Lagakos, M. A. Koch, C. Pettinelli, M. W. Myers, D. K. Booth, H. H. Balfour Jr, R. C. Reichman, J. A. Bartlett, M. S. Hirsch, et al. Zidovudine in asymptomatic human immunodeficiency virus infection: a controlled trial in persons with fewer than 500 cd4-positive cells per cubic mi...

  4. [1990]

    URLhttps://doi.org/10.1056/NEJM199004053221401. 32 A. Link selection in binary regression with the power prior When employing generalized linear models, the choice of link function can be crucial, as it controls the degree of non-linearity between the conditional mean and the linear predictor. In this section we illustrate how to employ hdbayesto compute ...

  5. [1995]

    URLhttps://doi.org/ 10.1080/01621459.1995.10476572. W. M. Landau.instantiate: Pre-Compiled ‘CmdStan’ Models in R Pack- ages,

  6. [1999]

    URLhttps://doi.org.fiss.iaha.eu/10.1111/ 1467-9868.00173. Y. Duan, K. Ye, and E. P. Smith. Evaluating water quality using power priors to incorporate historical information.Environmetrics: The Official Journal of the International Environmetrics Society, 17(1):95–106,

  7. [2000]

    URLhttps://doi.org/ 10.1214/ss/1009212673. J. G. Ibrahim, M.-H. Chen, Y. Gwon, and F. Chen. The power prior: theory and applications.Statistics in Medicine, 34(28):3724–3749,

  8. [2003]

    URLhttps://doi.org/10.1214/aos/1056562461. B. Neuenschwander, M. Branson, and D. J. Spiegelhalter. A note on the power prior.Statistics in Medicine, 28(28):3562–3566,

Show all 15 references
  1. [2006]

    URLhttps://doi.org/10.1002/env.752. M. Friedman. Piecewise exponential models for survival data with covariates. The Annals of Statistics, 10(1):101–113,

  2. [2009]

    R Core Team.R: A Language and Environment for Statistical Computing

    URLhttps: //doi.org/10.1002/sim.3722. R Core Team.R: A Language and Environment for Statistical Computing. R Foundation for Statistical Computing, Vienna, Austria,

  3. [2012]

    URLhttps://doi.org/10.1214/12-BA722. M. D. Hoffman, A. Gelman, et al. The No-U-Turn sampler: adaptively setting path lengths in Hamiltonian Monte Carlo.J. Mach. Learn. Res., 15 (1):1593–1623,

  4. [2014]

    URLhttps://doi.org/10.1111/biom.12242. M. Secrest and I. Gravestock.psborrow2: Bayesian Dynamic Borrowing Analysis and Simulation,

  5. [2015]

    URL https://doi.org/10.1002/sim.6728. R. E. Kass and A. E. Raftery. Bayes factors.Journal of the American Statistical Association, 90(430):773–795,

  6. [2017]

    URL https://doi.org/10.18637/jss.v076.i01. L. M. Carvalho and J. G. Ibrahim. On the normalized power prior.Statistics in Medicine, 40(24):5251–5275,

  7. [2020]

    URLhttps://doi.org/10.18637/jss.v092. i10. B. P. Hobbs, D. J. Sargent, and B. P. Carlin. Commensurate priors for incorporating historical information in clinical trials using general and generalized linear models.Bayesian Analysis (Online), 7(3):639,

Pith tools

Reviewed August 6, 2026 · model on record in the stance chip above.