nbcanlink {VGAM}

Negative Binomial Canonical Link Function

Description

Computes the negative binomial canonical link transformation, including its inverse and the first two derivatives.

Usage

nbcanlink(theta, size = NULL, wrt.param = NULL, bvalue = NULL,
          inverse = FALSE, deriv = 0, short = TRUE, tag = FALSE)

Arguments

theta

Numeric or character. Typically the mean of a negative binomial distribution (NBD). See below for further details.

size, wrt.param

size contains the k matrix which must be of a conformable dimension as theta. Also, if deriv > 0 then wrt.param is either 1 or 2 (1 for with respect to the first parameter, and 2 for with respect to the second parameter (size)).

bvalue

Details at Links.

inverse, deriv, short, tag

Details at Links.

Details

The NBD canonical link is \log(\theta/(\theta + k)) where \theta is the NBD mean. The canonical link is used for theoretically relating the NBD to GLM class.

This link function was specifically written for negbinomial and negbinomial.size, and should not be used elsewhere (these VGAM family functions have code that specifically handles nbcanlink().)

Estimation with the NB canonical link has a somewhat interesting history. If we take the problem as beginning with the admission of McCullagh and Nelder (1983; first edition, p.195) [see also McCullagh and Nelder (1989, p.374)] that the NB is little used in applications and has a “problematical” canonical link then it appears only one other publicized attempt was made to solve the problem seriously. This was Hilbe, who produced a defective solution. However, Miranda and Yee (2023) solve this four-decade old problem using total derivatives and it is implemented by using nbcanlink with negbinomial. Note that early versions of VGAM had a defective solution.

Value

For deriv = 0, the above equation when inverse = FALSE, and if inverse = TRUE then kmatrix / expm1(-theta) where theta is really eta. For deriv = 1, then the function returns d eta / d theta as a function of theta if inverse = FALSE, else if inverse = TRUE then it returns the reciprocal.

Note

While theoretically nice, this function is not recommended in general since its value is always negative (linear predictors ought to be unbounded in general). A loglink link for argument lmu is recommended instead.

Numerical instability may occur when theta is close to 0 or 1. Values of theta which are less than or equal to 0 can be replaced by bvalue before computing the link function value. See Links.

Author(s)

Victor Miranda and Thomas W. Yee.

References

Hilbe, J. M. (2011). Negative Binomial Regression, 2nd Edition. Cambridge: Cambridge University Press.

McCullagh, P. and Nelder, J. A. (1989). Generalized Linear Models, 2nd ed. London: Chapman & Hall.

Miranda-Soberanis, V. F. and Yee, T. W. (2023). Two-parameter link functions, with applications to negative binomial, Weibull and quantile regression. Computational Statistics, 38, 1463–1485.

Yee, T. W. (2014). Reduced-rank vector generalized linear models with two linear predictors. Computational Statistics and Data Analysis, 71, 889–902.

See Also

negbinomial, negbinomial.size.

Examples

nbcanlink("mu", short = FALSE)

mymu <- 1:10  # Test some basic operations:
kmatrix <- cbind(runif(length(mymu)))
eta1 <- nbcanlink(mymu, size = kmatrix)
ans2 <- nbcanlink(eta1, size = kmatrix, inverse = TRUE)
max(abs(ans2 - mymu))  # Should be 0

## Not run:  mymu <- seq(0.5, 10, length = 101)
kmatrix <- matrix(10, length(mymu), 1)
plot(nbcanlink(mymu, size = kmatrix) ~ mymu, las = 1,
     type = "l", col = "blue", xlab = expression({mu}))

## End(Not run)

# Estimate the parameters from some simulated data
ndata <- data.frame(x2 = runif(nn <- 500))
ndata <- transform(ndata, eta1 = -1 - 1 * x2,  # eta1 < 0
                          size1 = exp(1),
                          size2 = exp(2))
ndata <- transform(ndata,
            mu1 = nbcanlink(eta1, size = size1, inverse = TRUE),
            mu2 = nbcanlink(eta1, size = size2, inverse = TRUE))
ndata <- transform(ndata, y1 = rnbinom(nn, mu = mu1, size1),
                          y2 = rnbinom(nn, mu = mu2, size2))
summary(ndata)

nbcfit <-
  vglm(cbind(y1, y2) ~ x2,  #  crit = "c",
       negbinomial(lmu = "nbcanlink"),
       data = ndata, trace = TRUE)
coef(nbcfit, matrix = TRUE)
summary(nbcfit)