sampling_multinom {multinomineq}

Posterior Sampling for Inequality-Constrained Multinomial Models

Description

Uses Gibbs sampling to draw posterior samples for binomial and multinomial models with linear inequality-constraints.

Usage

sampling_multinom(
  k,
  options,
  A,
  b,
  V,
  prior = rep(1, sum(options)),
  M = 5000,
  start,
  burnin = 10,
  progress = TRUE,
  cpu = 1
)

sampling_binom(
  k,
  n,
  A,
  b,
  V,
  map = 1:ncol(A),
  prior = c(1, 1),
  M = 5000,
  start,
  burnin = 10,
  progress = TRUE,
  cpu = 1
)

Arguments

k

the number of choices for each alternative ordered by item type (e.g. c(a1,a2,a3, b1,b2) for a ternary and a binary item type). The length of k must be equal to the sum of options. The default k=0 is equivalent to sampling from the prior.

options

number of observable categories/probabilities for each item type/multinomial distribution, e.g., c(3,2) for a ternary and binary item.

A

a matrix with one row for each linear inequality constraint and one column for each of the free parameters. The parameter space is defined as all probabilities x that fulfill the order constraints A*x <= b.

b

a vector of the same length as the number of rows of A.

V

a matrix of vertices (one per row) that define the polytope of admissible parameters as the convex hull over these points (if provided, A and b are ignored). Similar as for A, columns of V omit the last value for each multinomial condition (e.g., a1,a2,a3,b1,b2 becomes a1,a2,b1). Note that this method is comparatively slow since it solves linear-programming problems to test whether a point is inside a polytope (Fukuda, 2004) or to run the Gibbs sampler.

prior

the prior parameters of the Dirichlet-shape parameters. Must have the same length as k.

M

number of posterior samples

start

only relevant if steps is defined or cmin>0: a vector with starting values in the interior of the polytope. If missing, an approximate maximum-likelihood estimate is used.

burnin

number of burnin samples that are discarded. Can be chosen to be small if the maxmimum-a-posteriori estimate is used as the (default) starting value.

progress

whether a progress bar should be shown (if cpu=1).

cpu

either the number of CPUs using separate MCMC chains in parallel, or a parallel cluster (e.g., cl <- parallel::makeCluster(3)). All arguments of the function call are passed directly to each core, and thus the total number of samples is M*number_cpu.

n

the number of choices per item type. If k=n=0, Bayesian inference is relies on the prior distribution only.

map

optional: numeric vector of the same length as k with integers mapping the frequencies k to the free parameters/columns of A/V, thereby allowing for equality constraints (e.g., map=c(1,1,2,2)). Reversed probabilities 1-p are coded by negative integers. Guessing probabilities of .50 are encoded by zeros. The default assumes different parameters for each item type: map=1:ncol(A)

Details

Draws posterior samples for binomial/multinomial random utility models that assume a mixture over predefined preference orders/vertices that jointly define a convex polytope via the set of inequalities A * x < b or as the convex hull of a set of vertices V.

Value

an mcmc matrix (or an mcmc.list if cpu>1) with posterior samples of the binomial/multinomial probability parameters. See mcmc) .

References

Myung, J. I., Karabatsos, G., & Iverson, G. J. (2005). A Bayesian approach to testing decision making axioms. Journal of Mathematical Psychology, 49, 205-225. doi:10.1016/j.jmp.2005.02.004

See Also

count_multinom, ml_multinom

Examples

############### binomial ##########################
A <- matrix(
  c(
    1, 0, 0, # x1 < .50
    1, 1, 1, # x1+x2+x3 < 1
    0, 2, 2, # 2*x2+2*x3 < 1
    0, -1, 0, # x2 > .2
    0, 0, 1
  ), # x3 < .1
  ncol = 3, byrow = TRUE
)
b <- c(.5, 1, 1, -.2, .1)
samp <- sampling_binom(c(5, 12, 7), c(20, 20, 20), A, b)
head(samp)
plot(samp)


############### multinomial ##########################
# binary and ternary choice:
#           (a1,a2   b1,b2,b3)
k <- c(15, 9, 5, 2, 17)
options <- c(2, 3)

# columns:   (a1,  b1,b2)
A <- matrix(
  c(
    1, 0, 0, # a1 < .20
    0, 2, 1, # 2*b1+b2 < 1
    0, -1, 0, # b1 > .2
    0, 0, 1
  ), # b2 < .4
  ncol = 3, byrow = TRUE
)
b <- c(.2, 1, -.2, .4)
samp <- sampling_multinom(k, options, A, b)
head(samp)
plot(samp)