Fits binomial-response GLMs using the bias-reduction method developed in Firth (1993) for the removal of the leading (\(\mathop{\rm O}(n^{-1})\)) term from the asymptotic expansion of the bias of the maximum likelihood estimator. Fitting is performed using pseudo-data representations, as described in Kosmidis (2007, Chapter 5). For estimation in binomial-response GLMs, the bias-reduction method is an improvement over traditional maximum likelihood because:
the bias-reduced estimator is second-order unbiased and has smaller variance than the maximum likelihood estimator and
the resulting estimates and their corresponding standard errors are always finite while the maximum likelihood estimates can be infinite (in situations where complete or quasi separation occurs).
brglm(formula, family = binomial, data, weights, subset, na.action,
start = NULL, etastart, mustart, offset,
control.glm = glm.control1(...), model = TRUE, method = "brglm.fit",
pl = FALSE, x = FALSE, y = TRUE, contrasts = NULL,
control.brglm = brglm.control(...), ...)brglm.fit(x, y, weights = rep(1, nobs), start = NULL, etastart = NULL,
mustart = NULL, offset = rep(0, nobs), family = binomial(),
control = glm.control(), control.brglm = brglm.control(),
intercept = TRUE, pl = FALSE)
as in glm.
as in glm. brglm currently
supports only the "binomial" family with links
"logit", "probit", "cloglog", "cauchit".
as in glm.
as in glm.
as in glm.
as in glm.
as in glm.
as in glm.
as in glm.
as in glm.
same as in glm. Only available to
brglm.fit.
as in glm.
as in glm.
the method to be used for fitting the model. The default
method is "brglm.fit", which uses either the modified-scores
approach to estimation or maximum penalized likelihood (see
the pl argument below). The standard
glm methods "glm.fit" for maximum likelihood and
"model.frame" for returning the model frame without any
fitting, are also accepted.
a logical value indicating whether the
model should be fitted using maximum penalized likelihood, where the
penalization is done using Jeffreys invariant prior, or using the
bias-reducing modified scores. It is only used when
method = "brglm.fit". The default value is FALSE (see also the
Details section).
as in glm.
as in glm.
as in glm.
a list of parameters for controlling the fitting
process when method = "brglm.fit". See documentation of
brglm.control for details.
further arguments passed to or from other methods.
brglm returns an object of class "brglm". A
"brglm" object inherits first from "glm" and then from
"lm" and is a list containing the following components:
as in glm.
as in glm.
as in glm.
as in glm.
as in glm.
as in glm.
as in glm.
as in glm.
as in glm.
as in glm.
as in glm (see Details).
as in glm.
as in glm.
as in glm.
as in glm.
as in glm.
as in glm.
as in glm.
as in glm.
as in glm.
the vector of the modified scores for the
parameters at the final iteration. If pl = TRUE they are the
derivatives of the penalized likelihood at the final iteration.
the Fisher information matrix evaluated at the
resulting estimates. Only available when method = "brglm.fit".
the diagonal elements of the hat matrix. Only available
when method = "brglm.fit"
the number of iterations that were required until
convergence. Only available when method = "brglm.fit".
a list with components ar and at which
contains the values of the additive modifications to the responses
(y) and to the binomial totals (prior.weights) at
the resulting estimates (see modifications for more
information). Only available when method = "brglm.fit".
as in glm.
as in glm.
as in glm.
as in glm.
as in glm.
as in glm.
as control in the result of
glm.
the control.brglm argument that was passed to
brglm. Only available when method = "brglm.fit".
the method used for fitting the model.
as in glm.
as in glm.
logical having the same value with the pl
argument passed to brglm. Only available when method =
"brglm.fit".
1. It is not advised to use methods associated with model comparison
(add1, drop1,
anova, etc.) on objects of class
"brglm". Model comparison when estimation is performed using
the modified scores or the penalized likelihood is an on-going
research topic and will be implemented as soon as it is concluded.
2. The use of Akaike's information criterion (AIC) for model selection
when method = "brglm.fit" is asymptotically valid, because
the log-likelihood derivatives dominate the modification (in terms
of asymptotic order).
brglm.fit is the workhorse function for fitting the model using
either the bias-reduction method or maximum penalized likelihood. If
method = "glm.fit", usual maximum likelihood is used via
glm.fit.
The main iteration of brglm.fit consists of the following
steps:
Calculate the diagonal components of the hat matrix (see
gethats and hatvalues).
Obtain the pseudo-data representation at the current value of the
parameters (see modifications for more information).
Fit a local GLM, using glm.fit on the pseudo data.
Adjust the quadratic weights to agree with the original binomial totals.
Iteration is repeated until either the iteration limit has been reached
or the sum of the absolute values of the modified scores is less than
some specified positive constant (see the br.maxit and
br.epsilon arguments in brglm.control).
The default value (FALSE) of pl, when method = "brglm.fit",
results in estimates that are free of any \(\mathop{\rm
O}(n^{-1})\) terms in the asymptotic expansion of their bias. When
pl = TRUE bias-reduction is again achieved but generally not at
such order of magnitude. In the case of logistic regression the value of
pl is irrelevant since maximum penalized likelihood and the
modified-scores approach coincide for natural exponential families (see
Firth, 1993).
For other language related details see the details section in
glm.
Kosmidis, I. and Firth, D. (2020). Jeffreys-prior penalty, finiteness and shrinkage in binomial-response generalized linear models. Biometrika. https://doi.org/10.1093/biomet/asaa052
Bull, S. B., Lewinger, J. B. and Lee, S. S. F. (2007). Confidence intervals for multinomial logistic regression in sparse data. Statistics in Medicine 26, 903--918.
Firth, D. (1992) Bias reduction, the Jeffreys prior and GLIM. In Advances in GLIM and statistical modelling: Proceedings of the GLIM 92 conference, Munich, Eds. L.~Fahrmeir, B.~Francis, R.~Gilchrist and G.Tutz, pp. 91--100. New York: Springer.
Firth, D. (1992) Generalized linear models and Jeffreys priors: An iterative generalized least-squares approach. In Computational Statistics I, Eds. Y. Dodge and J. Whittaker. Heidelberg: Physica-Verlag.
Firth, D. (1993). Bias reduction of maximum likelihood estimates. Biometrika 80, 27--38.
Heinze, G. and Schemper, M. (2002). A solution to the problem of separation in logistic regression. Statistics in Medicine 21, 2409--2419.
Kosmidis, I. (2007). Bias reduction in exponential family nonlinear models. PhD Thesis, Department of Statistics, University of Warwick.
Kosmidis, I. and Firth, D. (2009). Bias reduction in exponential family nonlinear models. Biometrika 96, 793--804.
# NOT RUN {
## Begin Example
data(lizards)
# Fit the GLM using maximum likelihood
lizards.glm <- brglm(cbind(grahami, opalinus) ~ height + diameter +
light + time, family = binomial(logit), data=lizards,
method = "glm.fit")
# Now the bias-reduced fit:
lizards.brglm <- brglm(cbind(grahami, opalinus) ~ height + diameter +
light + time, family = binomial(logit), data=lizards,
method = "brglm.fit")
lizards.glm
lizards.brglm
# Other links
update(lizards.brglm, family = binomial(probit))
update(lizards.brglm, family = binomial(cloglog))
update(lizards.brglm, family = binomial(cauchit))
# Using penalized maximum likelihood
update(lizards.brglm, family = binomial(probit), pl = TRUE)
update(lizards.brglm, family = binomial(cloglog), pl = TRUE)
update(lizards.brglm, family = binomial(cauchit), pl = TRUE)
# }
Run the code above in your browser using DataLab