Sample a binomial regression rate

Sample a binomial regression rate

sample_binom_reg(
  p,
  k,
  n,
  mean,
  precision,
  method = c("slice", "normal", "uniform", "quadratic taylor", "mv quadratic taylor",
    "mv ind quadratic taylor", "mv beta", "mv truncated exponential"),
  ...,
  width = 1,
  nexpand = 10,
  ncontract = 100,
  acceptance = c("MH", "LL only", "regardless")
)

Arguments

p

the previous iteration of the logit-probability

k

the realized value from the binomial distribution

n

the number of trials

mean

the prior mean. If this is negative infinity, precision must not be a matrix, and k must be 0.

precision

the prior precision

method

The method to use to propose new values.

  • "normal" proposes normal deviates from the current value.

  • "uniform" proposes uniform deviates.

  • "gamma" (only for Poisson) does moment-matching on the mean and variance to propose a conjugate gamma proposal.

  • "mv gamma" (only for Poisson) does the same, but ignores off-diagonal elements of the precision matrix in the proposal (for speed's sake; this may result in low acceptance rates when elements are highly correlated).

  • "mv beta" (only for binomial, and still somewhat experimental) does moment matching on a mean and variance approximated by the log-normal distribution (instead of logit-normal, for which there is no closed form) to propose a conjugate beta proposal; here again we ignore off-diagonal elements of the precision matrix for the proposal.

  • "quadratic taylor" proposes using a second-order Taylor approximation of the log-density at the current value, which amounts to a normal proposal with mean (usually) not equal to the current value.

  • "mv quadratic taylor" does the same, but uses a multivariate normal approximation.

  • "mv ind quadratic taylor" proposes using a similar Taylor approximation, but approximates the log-density around the mean, instead of the current value. Furthermore, like "mv beta" and "mv gamma", it ignores off-diagonal elements of the precision matrix for the proposal, which again might yield low acceptance rates when elements are highly correlated; in particular, it is not recommended for sample_multinom_reg, whose data are a priori correlated (unless you increase the variance a lot, in which case it works okay). Both of these simplifications yield a significant speed boost.

  • "mv truncated exponential" proposes using a first-order Taylor approximation of the log-density at the current value, which amounts to a truncated exponential proposal.

Note that "mv gamma", "mv beta" and "mv [ind ]quadratic taylor" accept or reject an entire row at a time.

...

Other arguments (not used)

width

For "normal" proposals, the standard deviation(s) of proposals. For "uniform" proposals, the half-width of the uniform proposal interval. For "slice", the width of each expansion (to the right and left each). For "gamma", "beta", and "mv ind quadratic taylor" a scaling factor to increase the variance of the proposal (it's squared, so that width is on the sd scale).

nexpand

The maximum number of expansions (to the right and left each)

ncontract

The maximum number of contractions. If this is exceeded, the original value is returned

acceptance

What should be the criteria for acceptance? "MH" indicates the usual Metropolis-Hastings update. "LL only" ignores the proposal densities but considers the log-likelihoods. "regardless" accepts no matter what. This is useful for testing, or for when the method is a gamma, beta, or quadratic approximation, which can be hard to accept if the initial starting point is low-density.

Details

This function samples p conditional on k, n, mean, and precision, where k ~ Binom(n, expit(p)) and p ~ N(mean, precision).

In the case that n is zero, sampling is ignored in favor of a normal draw. In the special case when an entire (multivariate) row of n is zero, the entire row is simultaneously drawn (in R); when only some elements are zero, they're drawn univariately (in C++). Note that only the latter is used when proposal indicates multivariate MH.

This is vectorized over p, k, n, and mean. If precision is a matrix, p is assumed to be multivariately distributed, and different internals are used.

The internals are defined in C++.