Distribution functions of the (inflated) ordered response

Probability mass (diord), cumulative probability (piord), quantile (qiord), and random generation (riord) for the ordered outcome of the models in this package, given the outcome linear predictor eta, the cutpoints tau, and – for inflated models – the inflation linear predictor(s) a, the inflated category k, and the error correlation rho; or taken from a fitted model (object, optionally at newdata). Categories are indexed 0, ..., J - 1 (the order of the fitted levels).

Usage

diord(
  x,
  eta = NULL,
  tau = NULL,
  a = NULL,
  k = NULL,
  link = c("probit", "logit"),
  rho = 0,
  object = NULL,
  newdata = NULL,
  log = FALSE
)

piord(
  q,
  eta = NULL,
  tau = NULL,
  a = NULL,
  k = NULL,
  link = c("probit", "logit"),
  rho = 0,
  object = NULL,
  newdata = NULL,
  lower.tail = TRUE,
  log.p = FALSE
)

qiord(
  p,
  eta = NULL,
  tau = NULL,
  a = NULL,
  k = NULL,
  link = c("probit", "logit"),
  rho = 0,
  object = NULL,
  newdata = NULL,
  lower.tail = TRUE,
  log.p = FALSE
)

riord(
  n,
  eta = NULL,
  tau = NULL,
  a = NULL,
  k = NULL,
  link = c("probit", "logit"),
  rho = 0,
  object = NULL,
  newdata = NULL
)

Arguments

x, q: Vector of categories (0-based indices, or level labels when object is given).
eta: Outcome linear predictor \(x'\beta\) (numeric vector).
tau: Cutpoints, strictly increasing (length J - 1).
a: Inflation linear predictor \(z'\gamma\): a vector (common split) or an n x (J - 1) matrix (category-specific split); NULL for the plain ordered model.
k: The inflated category, a 0-based index (or a level label when object is given); NULL for the plain ordered model.
link: "probit" or "logit".
rho: Error correlation (probit only; a vector with one entry per split equation under a category-specific split).
object: Optional fitted "iord" object from which eta, tau, a, k, link, and rho are taken (fixed-effects fits only).
newdata: Optional data frame of covariate values for object (default: the estimation data).
log, log.p, lower.tail: As in the base distribution functions.
p: Vector of probabilities.
n: Number of draws (one per element of eta; a scalar n with a single eta gives n draws).

Value

diord: probabilities P(y = x); piord: P(y <= q); qiord: categories (0-based, or labels when object is given); riord: drawn categories. For object-based calls the length is the number of rows of newdata (or of the estimation data), recycling x/q/p.

Details

With a = NULL the distribution is the plain ordered probit or logit, \(P(y = j) = F(\tau_j - \eta) - F(\tau_{j-1} - \eta)\). With a vector a it is the inflated model of iop() / iol(), \(P(y = j) = F(a)\,\pi_j + 1\{j = k\}[1 - F(a)]\) (with the bivariate normal rectangle probabilities when rho != 0), and with an n x (J - 1) matrix a it is the category-specific split model (split = "category"), one column per non-inflated category in order. qiord() returns the smallest category whose cumulative probability reaches p; riord() draws categories by inversion. The functions are vectorized over eta (recycling x, q, p against it); riop() remains the data-generating simulator that also draws covariates.

Examples

## plain ordered probit, one observation
diord(0:2, eta = 0.3, tau = c(-0.5, 0.8))
#> [1] 0.2118554 0.4796071 0.3085375
piord(1, eta = 0.3, tau = c(-0.5, 0.8))
#> [1] 0.6914625
## zero-inflated ordered probit: a split probability of F(0.4) = 0.66
diord(0:2, eta = 0.3, tau = c(-0.5, 0.8), a = 0.4, k = 0)
#> [1] 0.4834329 0.3143449 0.2022222
diord(0:2, eta = 0.3, tau = c(-0.5, 0.8), a = 0.4, k = 0, rho = -0.5)
#> [1] 0.5324408 0.3358691 0.1316901
qiord(c(0.1, 0.5, 0.9), eta = 0.3, tau = c(-0.5, 0.8), a = 0.4, k = 0)
#> [1] 0 1 2
table(riord(1000, eta = 0.3, tau = c(-0.5, 0.8), a = 0.4, k = 0))
#> 
#>   0   1   2 
#> 501 285 214 
## from a fitted model
data(bp)
m <- iop(violence ~ loggdppc + parliament + disaster | loggdppc + parliament, data = bp,
         inflate = "bottom")
#> The inflation equation contains no covariate that is excluded from the outcome equation; the split is then identified by functional form alone. An exclusion restriction is advisable.
head(diord("none", object = m))                      # P(y = none) for each observation
#> [1] 0.9633038 0.9461806 0.9729054 0.9808280 0.9567877 0.9755843
head(piord("repression", object = m))                # P(y <= repression)
#> [1] 0.9876523 0.9709297 0.9908706 0.9945956 0.9766412 0.9917695
qiord(0.5, object = m, newdata = bp[1:5, ])          # median category at five profiles
#> [1] none none none none none
#> Levels: none < repression < civil war