brmultinom: Bias reduction for multinomial response models using the Poisson trick.

Description

brmultinom is a wrapper of brglmFit that fits multinomial regression models using implicit and explicit bias reduction methods. See Kosmidis & Firth (2011) for details.

Usage

brmultinom(formula, data, weights, subset, na.action, contrasts = NULL,
  ref = 1, model = TRUE, x = TRUE, control = list(...), ...)

Arguments

formula

a formula expression as for regression models, of the form response ~ predictors. The response should be a factor or a matrix with K columns, which will be interpreted as counts for each of K classes. A log-linear model is fitted, with coefficients zero for the first class. An offset can be included: it should be a numeric matrix with K columns if the response is either a matrix with K columns or a factor with K >= 2 classes, or a numeric vector for a response factor with 2 levels. See the documentation of formula() for other details.

data

an optional data frame in which to interpret the variables occurring in formula.

weights

optional case weights in fitting.

subset

expression saying which subset of the rows of the data should be used in the fit. All observations are included by default.

na.action

a function to filter missing data.

contrasts

a list of contrasts to be used for some or all of the factors appearing as variables in the model formula.

ref

the reference category to use for multinomial regression. Either an integer, in which case levels(response)[ref] is used as a baseline, or a character string. Default is 1.

model

logical. If true, the model frame is saved as component model of the returned object.

should the model matrix be included with in the result (default is TRUE).

control

a list of parameters for controlling the fitting process. See brglmControl for details.

...

arguments to be used to form the default 'control' argument if it is not supplied directly.

Details

The models brmultinom handles are also known as baseline-category logit models (see, Agresti, 2002, Section 7.1), because they model the log-odds of every category against a baseline category. The user can control which baseline (or reference) category is used via the ref. By default brmultinom uses the first category as reference.

The maximum likelihood estimates for the parameters of baseline-category logit models have infinite components with positive probability, which can result in problems in their estimation and the use of inferential procedures (e.g. Wad tests). Albert and Andreson (1984) have categorised the possible data patterns for such models into the exclusive and exhaustive categories of complete separation, quasi-complete separation and overlap, and showed that infinite maximum likelihood estimates result when complete or quasi-complete separation occurs.

The adjusted score approach to bias reduction that brmultinom implements (type = "AS_mean") is an alternative to maximum likelihood that results in estimates with smaller asymptotic bias that are also *always* finite, even in cases of complete or quasi-complete separation.

brmultinom is a wrapper of brglmFit that fits multinomial logit regression models through the 'Poisson trick' (see, for example, Palmgren, 1981; Kosmidis & Firth, 2011).

The implementation relies on the construction of an 'extended' model matrix for the log-linear model and constraints on the sums of the Poisson means. Specifically, a log-linear model is fitted on a Kronecker product (https://en.wikipedia.org/wiki/Kronecker_product) of the original model matrix X implied by the formula, augmented by nrow(X) dummy variables.

The extended model matrix is sparse, and the Matrix package is used for its effective storage.

While brmultinom can be used for analyses using multinomial regression models, the current implementation is more of a 'proof of concept' and is not expected to scale well with either of nrow(X), ncol(X) or the number of levels in the categorical response.

References

Kosmidis I, Kenne Pagui EC, Sartori N (2019). Mean and median bias reduction in generalized linear models. *arXiv e-prints*, arXiv:1804.04085. To appear in Statistics and Computing, <URL: https://arxiv.org/abs/1804.04085>.

Agresti A. (2002). *Categorical data analysis* (2nd edition). Wiley Series in Probability and Statistics. Wiley.

Albert A. and Anderson J. A. (1984). On the Existence of Maximum Likelihood Estimates in Logistic Regression Models. *Biometrika*, **71** 1--10.

Kosmidis I. and Firth D. (2011). Multinomial logit bias reduction via the Poisson log-linear model. *Biometrika*, **98**, 755-759.

Palmgren, J. (1981). The Fisher Information Matrix for Log Linear Models Arguing Conditionally on Observed Explanatory Variables. *Biometrika*, **68**, 563-566.

Examples

Run this code

# NOT RUN {
data("housing", package = "MASS")

# Maximum likelihood using nnet::multinom
houseML1nnet <- nnet::multinom(Sat ~ Infl + Type + Cont, weights = Freq,
                               data = housing)
# Maximum likelihood using brmultinom with baseline category 'Low'
houseML1 <- brmultinom(Sat ~ Infl + Type + Cont, weights = Freq,
                       data = housing, type = "ML", ref = 1)
# The estimates are numerically the same as houseML0
all.equal(coef(houseML1nnet), coef(houseML1), tolerance = 1e-04)

# Maximum likelihood using brmultinom with 'High' as baseline
houseML3 <- brmultinom(Sat ~ Infl + Type + Cont, weights = Freq,
                      data = housing, type = "ML", ref = 3)
# The fitted values are the same as houseML1
all.equal(fitted(houseML3), fitted(houseML1), tolerance = 1e-10)

# Bias reduction
houseBR3 <- update(houseML3, type = "AS_mean")
# Bias correction
houseBC3 <- update(houseML3, type = "correction")


# }