ml_binom: Maximum-likelihood Estimate

View source: R/ml_estimates.R

ml_binomR Documentation

Maximum-likelihood Estimate

Description

Get ML estimate for product-binomial/multinomial model with linear inequality constraints.

Usage

ml_binom(k, n, A, b, map, strategy, n.fit = 3, start, progress = FALSE, ...)

ml_multinom(k, options, A, b, V, n.fit = 3, start, progress = FALSE, ...)

Arguments

k

vector of observed response frequencies.

n

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

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.

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)

strategy

a list that defines the predictions of a strategy, seestrategy_multiattribute.

n.fit

number of calls to constrOptim.

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.

progress

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

...

further arguments passed to the function constrOptim. To ensure high accuracy, the number of maximum iterations should be sufficiently large (e.g., by setting control = list(maxit = 1e6, reltol=.Machine$double.eps^.6), outer.iterations = 1000.

options

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

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.

Details

First, it is checked whether the unconstrained maximum-likelihood estimator (e.g., for the binomial: k/n) is inside the constrained parameter space. Only if this is not the case, nonlinear optimization with convex linear-inequality constrained is used to estimate (A) the probability parameters θ for the Ab-representation or (B) the mixture weights α for the V-representation.

Value

the list returned by the optimizer constrOptim, including the input arguments (e.g., k, options, A, V, etc.).

  • If the Ab-representation was used, par provides the ML estimate for the probability vector θ.

  • If the V-representation was used, par provides the estimates for the (usually not identifiable) mixture weights α that define the convex hull of the vertices in V, while p provides the ML estimates for the probability parameters. Because the weights must sum to one, the α-parameter for the last row of the matrix V is dropped. If the unconstrained ML estimate is inside the convex hull, the mixture weights α are not estimated and replaced by missings (NA).

Examples

# predicted linear order: p1 < p2 < p3 < .50
# (cf. WADDprob in ?strategy_multiattribute)
A <- matrix(
  c(
    1, -1, 0,
    0, 1, -1,
    0, 0, 1
  ),
  ncol = 3, byrow = TRUE
)
b <- c(0, 0, .50)
ml_binom(k = c(4, 1, 23), n = 40, A, b)[1:2]
ml_multinom(
  k = c(4, 36, 1, 39, 23, 17),
  options = c(2, 2, 2), A, b
)[1:2]


# probabilistic strategy:  A,A,A,B  [e1<e2<e3<e4<.50]
strat <- list(
  pattern = c(-1, -2, -3, 4),
  c = .5, ordered = TRUE, prior = c(1, 1)
)
ml_binom(c(7, 3, 1, 19), 20, strategy = strat)[1:2]


# vertex representation (one prediction per row)
V <- matrix(c(
  # strict weak orders
  0, 1, 0, 1, 0, 1, # a < b < c
  1, 0, 0, 1, 0, 1, # b < a < c
  0, 1, 0, 1, 1, 0, # a < c < b
  0, 1, 1, 0, 1, 0, # c < a < b
  1, 0, 1, 0, 1, 0, # c < b < a
  1, 0, 1, 0, 0, 1, # b < c < a

  0, 0, 0, 1, 0, 1, # a ~ b < c
  0, 1, 0, 0, 1, 0, # a ~ c < b
  1, 0, 1, 0, 0, 0, # c ~ b < a
  0, 1, 0, 1, 0, 0, # a < b ~ c
  1, 0, 0, 0, 0, 1, # b < a ~ c
  0, 0, 1, 0, 1, 0, # c < a ~ b

  0, 0, 0, 0, 0, 0 # a ~ b ~ c
), byrow = TRUE, ncol = 6)
ml_multinom(
  k = c(4, 1, 5, 1, 9, 0, 7, 2, 1), n.fit = 1,
  options = c(3, 3, 3), V = V
)

multinomineq documentation built on Nov. 22, 2022, 5:09 p.m.