f_gamma: Specify a Gamma sampling distribution

In mcmcsae: Markov Chain Monte Carlo Small Area Estimation

Specify a Gamma sampling distribution

Description

This function can be used in the family argument of create_sampler or generate_data to specify a Gamma sampling distribution.

Usage

f_gamma(
  link = "log",
  shape.vec = ~1,
  shape.prior = pr_gamma(0.1, 0.1),
  control = set_MH(type = "RWLN", scale = 0.2, adaptive = TRUE)
)

Arguments

link	the name of a link function. Currently the only allowed link function for the gamma distribution is "log".
shape.vec	optional formula specification of unequal shape parameter.
shape.prior	prior for gamma shape parameter. Supported prior distributions: pr_fixed with a default value of 1, pr_exp and pr_gamma. The current default is pr_gamma(shape=0.1, rate=0.1).
control	options for the Metropolis-Hastings algorithm employed in case the shape parameter is to be inferred. Function set_MH can be used to change the default options. The two choices of proposal distribution type supported are "RWLN" for a random walk proposal on the log-shape scale, and "gamma" for an approximating gamma proposal, found using an iterative algorithm. In the latter case, a Metropolis-Hastings accept-reject step is currently omitted, so the sampling algorithm is an approximate one, though often quite accurate and efficient.

Value

A family object.

References

J.W. Miller (2019). Fast and Accurate Approximation of the Full Conditional for Gamma Shape Parameters. Journal of Computational and Graphical Statistics 28(2), 476-480.

Examples

## Not run: 
n <- 3000
m <- 25
dat <- data.frame(
  x = rnorm(n),
  g = factor(sample(1:m, n, replace=TRUE), levels=1:m)
)
v <- rnorm(m, sd=0.6)
alpha <- 1
mu <- exp(with(dat, 1 - 0.5*x + v[g]))
dat$y <- rgamma(n, shape=alpha, rate=alpha/mu)

sampler <- create_sampler(
  y ~ reg(~ x, name="beta") +       # fixed effects
      gen(factor = ~ g, name="v"),  # random intercepts
  data=dat, family="gamma"
)
sim <- MCMCsim(sampler, store.all=TRUE)
compute_DIC(sim)
waic(sim)
summary(sim)

bayesplot::mcmc_recover_intervals(as.array(sim$gamma_shape_), alpha)
bayesplot::mcmc_recover_intervals(as.array(sim$beta), c(1, -0.5))
bayesplot::mcmc_recover_intervals(as.array(sim$v_sigma), 0.6)
yrep <- predict(sim, iters=sample(1:1000, 10))
bayesplot::pp_check(dat$y, as.matrix(yrep), bayesplot::ppc_dens_overlay) + ggplot2::xlim(0, 25)

## End(Not run)

mcmcsae documentation built on July 21, 2026, 1:06 a.m.