Skip to content

cran/glmbayes

BranchesTags

Folders and files

NameName
Last commit message
Last commit date

Latest commit

 

History

3 Commits
R
R
 
 
build
build
 
 
data
data
 
 
demo
demo
 
 
inst
inst
 
 
man
man
 
 
src
src
 
 
tests
tests
 
 
tools
tools
 
 
vignettes
vignettes
 
 
DESCRIPTION
DESCRIPTION
 
 
MD5
MD5
 
 
NAMESPACE
NAMESPACE
 
 
NEWS.md
NEWS.md
 
 
README.md
README.md
 
 
cleanup
cleanup
 
 
configure
configure
 
 
configure.win
configure.win
 
 

Repository files navigation

glmbayes

CRAN status CRAN downloads Monthly downloads GitHub release (latest by date) License: GPL-3 GitHub Workflow Status

glmbayes provides independent and identically distributed (iid) samples for Bayesian Generalized Linear Models (GLMs). Its primary interface, glmb(), serves as a Bayesian analogue to R's glm() function, supporting Gaussian, Poisson, Binomial, and Gamma families under log-concave likelihoods. Sampling for most models is performed using accept-reject methods based on likelihood subgradients (Nygren and Nygren, 2006). For Gaussian models, the package also includes lmb(), a Bayesian counterpart to R's lm().

The package includes a rich set of supporting tools for prior specification, model diagnostics, and method functions that mirror those for lm() and glm(). Most functions are extensively documented, and a comprehensive set of vignettes are available to guide users through the package's capabilities.

This repository is 0.9.6 in development. The current CRAN release is version 0.9.5 (CRAN). The GitHub repository holds the source; R-Universe builds binaries from it. See NEWS.md for changes.

Installation

CRAN (release 0.9.5)

install.packages("glmbayes")

GitHub / R-Universe (install from both CRAN and R-Universe repositories if you want R-Universe binaries or faster mirrors):

install.packages("glmbayes",
                 repos = c("https://cloud.r-project.org",
                           "https://knygren.r-universe.dev"))

Prebuilt binaries from CRAN (0.9.5) and R-Universe are built without OpenCL GPU support. For the CRAN release, OpenCL requires installing from source on a system with OpenCL development files available. To set up GPU acceleration, follow

Chapter 16 — Large models: GPU acceleration using OpenCL https://knygren.r-universe.dev/articles/glmbayes/Chapter-16.html

Minimal Working Example

library(glmbayes)

# Dobson (1990), p. 93: Randomized Controlled Trial
counts <- c(18,17,15,20,10,20,25,13,12)
outcome <- gl(3,1,9)
treatment <- gl(3,3)
print(d.AD <- data.frame(treatment, outcome, counts))

## Classical glm
glm.D93 <- glm(counts ~ outcome + treatment,
               family = poisson())

## Bayesian glmb
# Step 1: Set up prior
ps <- Prior_Setup(counts ~ outcome + treatment, family = poisson())
mu <- ps$mu
V  <- ps$Sigma

# Step 2: Fit using glmb
glmb.D93 <- glmb(counts ~ outcome + treatment,
                 family = poisson(),
                 pfamily = dNormal(mu = mu, Sigma = V))

summary(glmb.D93)

Supported families, links, and pfamilies

As with glm(), models are defined by a formula for the linear predictor and a family() describing the likelihood and link. In addition, glmb() requires a pfamily object specifying the prior.

Priors on regression coefficients

The primary table below covers priors on the regression coefficients β. The standard prior for all families is dNormal. The conjugate priors dBeta and dGamma(Inv_Dispersion = FALSE) provide closed-form IID posterior draws for intercept-only models with an identity link.

Likelihood family Link functions Compatible pfamilies (coefficient priors)
Gaussian identity dNormal, dNormal_Gamma, dIndependent_Normal_Gamma
Poisson / Quasi-Poisson log dNormal
Poisson identity (intercept-only) dGamma(Inv_Dispersion = FALSE) — conjugate Gamma–Poisson rate prior
Binomial / Quasi-Binomial logit, probit, cloglog dNormal
Binomial identity (intercept-only) dBeta — conjugate Beta–Binomial probability prior
Gamma log dNormal
Gamma identity (intercept-only) dGamma(Inv_Dispersion = FALSE) — conjugate Gamma–Gamma rate prior

dNormal_Gamma and dIndependent_Normal_Gamma also model precision jointly with the coefficients; see the precision/dispersion table below.

Priors on precision / dispersion

dGamma(Inv_Dispersion = TRUE) (the default when Inv_Dispersion is omitted) places a Gamma prior on the inverse dispersion 1/φ with the regression coefficients β held fixed. This is the precision prior used in Gibbs sampling steps for dispersion estimation.

Likelihood family Link Compatible pfamilies (precision prior)
Gaussian identity dGamma — prior on 1/σ² (precision)
Gamma log dGamma — prior on 1/φ (shape / dispersion)

dNormal_Gamma and dIndependent_Normal_Gamma model β and precision jointly in a single conjugate step, avoiding the need for a separate Gibbs precision update.

Prior_Setup

For a default, data‑aligned prior using the same formula and family as glm(), call Prior_Setup(formula, family, data = ..., ...). The returned list includes default settings for the following:

  • mu, Sigma — Zellner‑style normal prior components for use with most priors
  • Additional Gaussian‑specific calibration components:
    • dispersion for use with the dNormal() prior (gaussian and Gamma families)
    • Sigma_0, shape and rate for use with the dNormal_Gamma() prior
    • shape_ING and rate for use with dIndependent_Normal_Gamma() prior
    • shape, rate_gamma and coefficients for use with the dGamma() precision prior
  • Conjugate prior calibration components (intercept-only models):
    • conj_beta (shape1, shape2, beta) for use with dBeta() (Binomial/identity)
    • conj_poisson (shape, rate, beta) for use with dGamma(Inv_Dispersion = FALSE) (Poisson/identity)

Optional arguments adjust prior weight, centering, and related settings (see the function help and vignette Chapter 04).

Typical Prior_Setup wiring

Assuming ps <- Prior_Setup(...):

  • Non‑Gaussian families (log/logit/probit/cloglog links):
    Use dNormal(mu = ps$mu, Sigma = ps$Sigma).
    (For Gamma GLMs, also supply dispersion from the fitted GLM or from ps; see example("glmb").)

  • Binomial — conjugate Beta prior (identity link, intercept-only):
    Use dBeta(shape1 = ps$conj_beta$shape1, shape2 = ps$conj_beta$shape2, beta = ps$conj_beta$beta).

  • Poisson — conjugate Gamma rate prior (identity link, intercept-only):
    Use dGamma(shape = ps$conj_poisson$shape, rate = ps$conj_poisson$rate, beta = ps$conj_poisson$beta, Inv_Dispersion = FALSE).

  • Gaussian — normal prior with known dispersion:
    Use dNormal(mu = ps$mu, Sigma = ps$Sigma, dispersion = ps$dispersion).

  • Gaussian — conjugate Normal–Gamma:
    Use dNormal_Gamma(mu = ps$mu, Sigma_0 = ps$Sigma_0, shape = ps$shape, rate = ps$rate).

  • Gaussian — independent Normal–Gamma:
    Use dIndependent_Normal_Gamma(mu = ps$mu, Sigma = ps$Sigma, shape = ps$shape_ING, rate = ps$rate).

  • Gaussian / Gamma — precision prior (coefficients fixed, for Gibbs):
    With rate_dg <- if (!is.null(ps$rate_gamma)) ps$rate_gamma else ps$rate, use
    dGamma(shape = ps$shape, rate = rate_dg, beta = ps$coefficients).

The default priors have limiting behaviors that produce estimates resembling classical estimates as priors get weak (see documentation and vignettes for details).

All supported models have log‑concave likelihoods, enabling efficient iid sampling via enveloping functions and subgradient‑based accept–reject algorithms, especially for models lacking standard iid samplers.

Examples and Demos

Use example() and demo() to explore built-in examples and demos for supported families and links:

## Bayesian linear regression
example("lmb")

## Bayesian generalized linear models
example("glmb")

## Beta-Binomial conjugacy: dBeta() prior; Bechdel test (requires bayesrules)
## See also: vignette("Chapter-02-S03", package = "glmbayes")
demo("Ex_12_BetaBinomial")

## Gamma-Poisson conjugacy: dGamma(Inv_Dispersion=FALSE); bike counts + heart
## transplant mortality (requires bayesrules; Appendix A requires LearnBayes)
## See also: vignette("Chapter-02-S04", package = "glmbayes")
demo("Ex_13_GammaPoisson")

## Predictions for fitted glmb objects (newdata, type, etc.)
example("predict.glmb")

## Deviance residuals and simulate() for posterior predictive checks (menarche)
example("residuals.glmb")

## Two-block Gibbs sampler compared with iid sampling (linear model)
example("rlmb")

## Default prior specification using Prior_Setup
example("Prior_Setup")

## Matrix-input GLM example with an informative prior
example("rglmb")

## Two-step Boston example: estimates and summarizes models with unknown
## dispersion using dGamma priors via rGamma_reg, rglmb, rlmb, glmb, and lmb
example("summary.rGamma_reg")

## High-dimensional Gaussian model (14 predictors) with GPU acceleration (requires OpenCL)
example("Boston_centered")

## High-dimensional binomial model (14 predictors) with GPU acceleration (requires OpenCL)
example("Cleveland")

## Hierarchical linear model (Rubin/Gelman 8-schools) via rlmb
demo("Ex_07_Schools")

## Hierarchical generalized linear model (Poisson BikeSharing) via rglmb
demo("Ex_09_BikeSharingPoisson")

## Detailed simulation pipeline for rNormalGLM models (JASA 2006; Vignette Chapter A05)
example("rNormalGLM_std")

## Detailed simulation pipeline for rIndepNormalGammaReg models (Vignette Chapter A07)
example("rIndepNormalGammaReg_std")

Methodology

For generalized linear models where well known sampling methods are unavailable, sampling follows the framework from Nygren and Nygren (2006), using likelihood subgradients to construct enveloping functions for the posterior distribution. When the posterior is approximately normal, the expected number of draws per acceptance is bounded as per that paper and as discussed in our vignettes. Dispersion can be sampled via rGamma_reg() (standalone) or jointly with coefficients via rNormalGamma_reg() and rindepNormalGamma_reg().

GPU Acceleration Using OpenCL

The implemented algorithms tend to have acceptable performance on CPUs up to around 10-14 dimensions. For larger models, the envelope construction is embarrassingly parallel. To accelerate envelope construction in such cases, the package provides optional GPU acceleration using OpenCL. This requires that users have GPU enabled machines and an OpenCL installation. These features are discussed in more detail in two of our vignettes.

Vignettes

The glmbayes package includes a comprehensive set of vignettes organized into five major parts. These vignettes guide users from introductory material through applied modeling, advanced topics, and the underlying simulation methods that support the package.

Part 1: An Introduction

Overview of the package, its design philosophy, single-parameter conjugate models, and the basic workflow for fitting Bayesian linear and generalized linear models.

Part 2: Bayesian regression models

These chapters focus on Bayesian linear regression (Gaussian family). Topics include lmb() fitting, Prior_Setup(), posterior predictive checks (bayesplot), deviance residuals and model summaries, bayestestR-style summaries, and the bridge to Bayesian GLMs in Part 3.

Part 3: Generalized Linear Models

This part presents Bayesian GLMs across the major likelihood families, including binomial, quasi-binomial, Poisson, quasi-Poisson, and Gamma models. It covers model specification, link functions, log-concavity, diagnostics, interpretation of posterior results, and tooling (bayesplot, bayestestR) for visualization and summaries.

Part 4: Advanced Topics

These chapters explore more complex modeling scenarios and computational strategies, such as informative priors, two-block Gibbs sampling, linear and generalized linear mixed-effects models, models with unknown dispersion parameters, and large-scale model fitting using GPU acceleration using OpenCL.

Part 5: Simulation Methods and Technical Implementation

This part documents the mathematical and algorithmic foundations of the package. Topics include estimation procedures, likelihood subgradient densities, envelope construction, accept-reject sampling, and technical reports on sampler design including implementation aspects for GPU acceleration using OpenCL.

Together, these vignettes form a comprehensive reference that supports users at all levels, from first-time Bayesian GLM users to researchers interested in the mathematical and computational details behind the samplers.

Feature Highlights

  • S3 interface mirroring the structure of base glm()
  • Posterior predictive checks via pp_check() from the 'bayesplot' package for fitted glmb objects
  • Accept-reject sampling for log-concave likelihoods
  • Samplers for both fixed and variable dispersion
  • Extensive vignettes to guide users through the package's capabilities
  • Modular prior setup function

Limitations

  • Non-log-concave likelihoods are not currently supported

Future Plans

  • R Mathlib (nmath) usage from C: Today the package vendors local copies of selected R Mathlib routines and headers in *.c sources. The plan is to switch to calling the same nmath functions that ship with R, via the supported linking/API path, so maintenance tracks base R instead of duplicating sources.
  • OpenCL / GPU code upstream: Routines currently living under the openclport and nmathopencl namespaces are slated to move into dedicated upstream packages. nmathopencl is already available on R-Universe; a CRAN release is targeted, after which glmbayes can depend on that package for a substantial share of OpenCL- and GPU-related functionality rather than carrying those implementations here.
  • Conjugate priors for intercept-only GLMs: Add pfamily specifications that supply conjugate priors for intercept-only glm()-style models (a single mean structure / scalar linear predictor), complementing the existing prior families for general designs.
  • bayestestR integration: Add methods or small wrappers so bayestestR summaries and diagnostics can be used with glmb / lmb fits in the same way as with other Bayesian modeling workflows.

Further performance and algorithm work:

  • Poisson speed (OpenCL and simulation): Precompute the log-factorial term log(y!) once per observation and reuse it in both OpenCL envelope construction and accept-reject simulation, since it depends only on the response, to reduce redundant lgamma evaluation and improve performance for large Poisson models.
  • Grid selection (simulation): Precompute cumulative PLSD and use inverse CDF sampling (e.g. binary search) to select the grid component per candidate instead of scanning PLSD, improving the simulation loop when many candidates are evaluated.

About

❗ This is a read-only mirror of the CRAN R package repository. glmbayes — Bayesian Generalized Linear Models (IID Samples). Homepage: https://CRAN.R-project.org/package=glmbayes, https://github.com/knygren/glmbayes, https://knygren.r-universe.dev/glmbayes Report bugs for this package: https://github.com/k ...

Resources

Readme
Activity
Custom properties

Stars

0 stars

Watchers

0 watching

Forks

0 forks
Report repository

Releases

3 tags

Packages

 
 
 

Contributors

Languages