confintvglm: Confidence Intervals for Parameters of VGLMs

Description

Computes confidence intervals (CIs) for one or more parameters in a fitted model. Currently the object must be a "vglm" object.

Usage

confintvglm(object, parm, level = 0.95, method = c("wald", "profile"),
            trace = NULL, …)

Arguments

object

A fitted model object.

parm, level, …

Same as confint.

method

Character. The default is the first method. Abbreviations are allowed. Currently "profile" is basically working; and it is likely to be more accurate especially for small samples, as it is based on a profile log likelihood, however it is computationally intensive.

trace

Logical. If TRUE then one can monitor the computation as it progresses (because it is expensive). The default is the orginal model's trace value (see vglm.control). Setting FALSE suppresses all intermediate output.

Value

Same as confint.

Details

The default for this methods function is based on confint.default and assumes asymptotic normality. In particular, the coef and vcov methods functions are used for vglm-class objects.

When method = "profile" the function profilevglm is called to do the profiling. The code is very heavily based on profile.glm which was originally written by D. M. Bates and W. N. Venables (For S in 1996) and subsequently corrected by B. D. Ripley. Sometimes the profiling method can give problems, for example, cumulative requires the \(M\) linear predictors not to intersect in the data cloud. Such numerical problems are less common when method = "wald", however, it is well-known that inference based on profile likelihoods is generally more accurate than Wald, especially when the sample size is small. The deviance (deviance(object)) is used if possible, else the difference 2 * (logLik(object) - ell) is computed, where ell are the values of the loglikelihood on a grid.

For Wald CIs and rrvglm-class objects, currently an error message is produced because I haven't gotten around to write the methods function; it's not too hard, but am too busy! An interim measure is to coerce the object into a "vglm" object, but then the confidence intervals will tend to be too narrow because the estimated constraint matrices are treated as known.

For Wald CIs and vgam-class objects, currently an error message is produced because the theory is undeveloped.

Examples

Run this code

# NOT RUN {
# Example 1: this is based on a glm example
counts <- c(18,17,15,20,10,20,25,13,12)
outcome <- gl(3, 1, 9); treatment <- gl(3, 3)
 glm.D93 <-  glm(counts ~ outcome + treatment, family = poisson())
vglm.D93 <- vglm(counts ~ outcome + treatment, family = poissonff)
confint(glm.D93) # needs MASS to be present on the system
confint.default(glm.D93)  # based on asymptotic normality
confint(vglm.D93)
confint(vglm.D93) - confint(glm.D93)    # Should be all 0s
confint(vglm.D93) - confint.default(glm.D93)  # based on asympt. normality

# Example 2: simulated negative binomial data with multiple responses
ndata <- data.frame(x2 = runif(nn <- 100))
ndata <- transform(ndata, y1 = rnbinom(nn, mu = exp(3+x2), size = exp(1)),
                          y2 = rnbinom(nn, mu = exp(2-x2), size = exp(0)))
fit1 <- vglm(cbind(y1, y2) ~ x2, negbinomial, data = ndata, trace = TRUE)
coef(fit1)
coef(fit1, matrix = TRUE)
confint(fit1)
confint(fit1, "x2:1")  #  This might be improved to "x2" some day...
# }
# NOT RUN {
confint(fit1, method = "profile")  # Computationally expensive
confint(fit1, "x2:1", method = "profile", trace = FALSE)
# }
# NOT RUN {
fit2 <- rrvglm(y1 ~ x2, negbinomial(zero = NULL), data = ndata)
confint(as(fit2, "vglm"))  # Too narrow (SEs are biased downwards)
# }

Run the code above in your browser using DataCamp Workspace