multilogitlink.RdComputes the multilogit transformation, including its inverse and the first two derivatives.
multilogitlink(theta, refLevel = "(Last)", M = NULL, whitespace = FALSE,
bvalue = NULL, inverse = FALSE, deriv = 0, all.derivs = FALSE,
short = TRUE, tag = FALSE)Numeric or character. See below for further details.
See multinomial.
See Links.
Logical. This is currently experimental only.
Details at Links.
The multilogitlink() link function is a generalization of the
logitlink link to \(M\) levels/classes. It forms the
basis of the multinomial logit model. It is sometimes
called the multi-logit link or the multinomial logit
link; some people use softmax too. When its inverse function
is computed it returns values which are positive and add to unity.
For multilogitlink with deriv = 0,
the multilogit of theta,
i.e.,
log(theta[, j]/theta[, M+1]) when inverse = FALSE,
and if inverse = TRUE then
exp(theta[, j])/(1+rowSums(exp(theta))).
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.
Here, all logarithms are natural logarithms, i.e., to base e.
McCullagh, P. and Nelder, J. A. (1989). Generalized Linear Models, 2nd ed. London: Chapman & Hall.
Numerical instability may occur when theta is
close to 1 or 0 (for multilogitlink).
One way of overcoming this is to use, e.g., bvalue.
Currently care.exp() is used to avoid NAs being
returned if the probability is too close to 1.
pneumo <- transform(pneumo, let = log(exposure.time))
fit <- vglm(cbind(normal, mild, severe) ~ let, # For illustration only!
multinomial, trace = TRUE, data = pneumo)
#> Iteration 1: deviance = 5.407271
#> Iteration 2: deviance = 5.34745
#> Iteration 3: deviance = 5.347382
#> Iteration 4: deviance = 5.347382
fitted(fit)
#> normal mild severe
#> 1 0.9927503 0.005875947 0.001373768
#> 2 0.9329702 0.043219077 0.023810688
#> 3 0.8488899 0.085745054 0.065365011
#> 4 0.7485338 0.128835331 0.122630879
#> 5 0.6393787 0.168725388 0.191895881
#> 6 0.5334715 0.201127232 0.265401245
#> 7 0.4313692 0.226188995 0.342441766
#> 8 0.3581471 0.239824757 0.402028109
predict(fit)
#> log(mu[,1]/mu[,3]) log(mu[,2]/mu[,3])
#> 1 6.5829217 1.45330985
#> 2 3.6682387 0.59614746
#> 3 2.5639424 0.27139129
#> 4 1.8089375 0.04935622
#> 5 1.2035440 -0.12868047
#> 6 0.6981629 -0.27730511
#> 7 0.2308628 -0.41473071
#> 8 -0.1155781 -0.51661353
multilogitlink(fitted(fit))
#> [,1] [,2]
#> 1 6.5829217 1.45330985
#> 2 3.6682387 0.59614746
#> 3 2.5639424 0.27139129
#> 4 1.8089375 0.04935622
#> 5 1.2035440 -0.12868047
#> 6 0.6981629 -0.27730511
#> 7 0.2308628 -0.41473071
#> 8 -0.1155781 -0.51661353
multilogitlink(fitted(fit)) - predict(fit) # Should be all 0s
#> [,1] [,2]
#> 1 0.000000e+00 2.220446e-16
#> 2 0.000000e+00 -1.110223e-16
#> 3 0.000000e+00 -5.551115e-17
#> 4 0.000000e+00 1.387779e-17
#> 5 0.000000e+00 -8.326673e-17
#> 6 0.000000e+00 0.000000e+00
#> 7 1.665335e-16 0.000000e+00
#> 8 -4.163336e-17 1.110223e-16
multilogitlink(predict(fit), inverse = TRUE) # rowSums() add to unity
#> [,1] [,2] [,3]
#> 1 0.9927503 0.005875947 0.001373768
#> 2 0.9329702 0.043219077 0.023810688
#> 3 0.8488899 0.085745054 0.065365011
#> 4 0.7485338 0.128835331 0.122630879
#> 5 0.6393787 0.168725388 0.191895881
#> 6 0.5334715 0.201127232 0.265401245
#> 7 0.4313692 0.226188995 0.342441766
#> 8 0.3581471 0.239824757 0.402028109
multilogitlink(predict(fit), inverse = TRUE, refLevel = 1)
#> [,1] [,2] [,3]
#> 1 0.001373768 0.9927503 0.005875947
#> 2 0.023810688 0.9329702 0.043219077
#> 3 0.065365011 0.8488899 0.085745054
#> 4 0.122630879 0.7485338 0.128835331
#> 5 0.191895881 0.6393787 0.168725388
#> 6 0.265401245 0.5334715 0.201127232
#> 7 0.342441766 0.4313692 0.226188995
#> 8 0.402028109 0.3581471 0.239824757
multilogitlink(predict(fit), inverse = TRUE) -
fitted(fit) # Should be all 0s
#> [,1] [,2] [,3]
#> 1 0 0 0
#> 2 0 0 0
#> 3 0 0 0
#> 4 0 0 0
#> 5 0 0 0
#> 6 0 0 0
#> 7 0 0 0
#> 8 0 0 0
multilogitlink(fitted(fit), deriv = 1)
#> normal mild severe
#> 1 138.943754 171.191239 728.926256
#> 2 15.990591 24.183102 43.022338
#> 3 7.795702 12.756267 16.368641
#> 4 5.312622 8.909735 9.294324
#> 5 4.337011 7.129762 6.448624
#> 6 4.018006 6.223741 5.129167
#> 7 4.076810 5.713387 4.440982
#> 8 4.350138 5.485197 4.159708
multilogitlink(fitted(fit), deriv = 2)
#> normal mild severe
#> 1 19025.449830 -28962.03409 -5.298736e+05
#> 2 221.420110 -534.27143 -1.762778e+03
#> 3 42.406154 -134.81708 -2.329056e+02
#> 4 14.029214 -58.92861 -6.519766e+01
#> 5 5.243333 -33.67970 -2.562486e+01
#> 6 1.080754 -23.15364 -1.234382e+01
#> 7 -2.281339 -17.87591 -6.214829e+00
#> 8 -5.368762 -15.65599 -3.390448e+00