CommonVGAMffArguments: Common VGAM Family Function Arguments

Description

Here is a description of some common and typical arguments found in many VGAM family functions, e.g., zero, lsigma, isigma, gsigma, eq.mean, nsimEI and parallel.

Usage

TypicalVGAMfamilyFunction(lsigma = "loglink",
                          isigma = NULL,
                          link.list = list("(Default)" = "identitylink",
                                           x2          = "loglink",
                                           x3          = "logofflink",
                                           x4          = "multilogitlink",
                                           x5          = "multilogitlink"),
                          earg.list = list("(Default)" = list(),
                                           x2          = list(),
                                           x3          = list(offset = -1),
                                           x4          = list(),
                                           x5          = list()),
                          gsigma = exp(-5:5),
                          eq.mean = FALSE,
                          parallel = TRUE,
                          ishrinkage = 0.95,
                          nointercept = NULL, imethod = 1,
                          type.fitted = c("mean", "quantiles", "Qlink",
                                          "pobs0", "pstr0", "onempstr0"),
                          percentiles = c(25, 50, 75),
                          probs.x = c(0.15, 0.85),
                          probs.y = c(0.25, 0.50, 0.75),
                          multiple.responses = FALSE, earg.link = FALSE,
                          whitespace = FALSE, bred = FALSE, lss = TRUE,
                          oim = FALSE, nsimEIM = 100, byrow.arg = FALSE,
                          nrfs = 1, zero = NULL)

Value

An object of class "vglmff" (see vglmff-class). The object is used by modelling functions such as vglm

and vgam.

Arguments

lsigma

Character. Link function applied to a parameter and not necessarily a mean. See Links for a selection of choices. If there is only one parameter then this argument is often called link.

link.list, earg.list: Some VGAM family functions (such as normal.vcm) implement models with potentially lots of parameter link functions. These two arguments allow many such links and extra arguments to be inputted more easily. One has something like link.list = list ("(Default)" = "identitylink", x2 = "loglink", x3 = "logofflink") and earg.list = list ("(Default)" = list(), x2 = list(), x3 = "list(offset = -1)"). Then any unnamed terms will have the default link with its corresponding extra argument. Note: the multilogitlink link is also possible, and if so, at least two instances of it are necessary. Then the last term is the baseline/reference group.
isigma: Optional initial values can often be inputted using an argument beginning with "i". For example, "isigma" and "ilocation", or just "init" if there is one parameter. A value of NULL means a value is computed internally, i.e., a self-starting VGAM family function. If a failure to converge occurs make use of these types of arguments.
gsigma: Grid-search initial values can be inputted using an argument beginning with "g", e.g., "gsigma", "gshape" and "gscale". If argument isigma is inputted then that has precedence over gsigma, etc. If the grid search is 2-dimensional then it is advisable not to make the vectors too long as a nested for loop may be used. Ditto for 3-dimensions etc. Sometimes a ".mux" is added as a suffix, e.g., gshape.mux; this means that the grid is created relatively and not absolutely, e.g., its values are multipled by some single initial estimate of the parameter in order to create the grid on an absolute scale.
eq.mean: Logical. Constrain all the means to be equal? This type of argument is simpler than parallel because only a single TRUE or FALSE can be assigned; not a formula. Thus if TRUE then it will be enforced over all variables.
parallel: A logical, or a simple formula specifying which terms have equal/unequal coefficients. The formula must be simple, i.e., additive with simple main effects terms. Interactions and nesting etc. are not handled. To handle complex formulas use the constraints argument (of vglm etc.); however, there is a lot more setting up involved and things will not be as convenient.

nsimEIM: Some VGAM family functions use simulation to obtain an approximate expected information matrix (EIM). For those that do, the nsimEIM argument specifies the number of random variates used per observation; the mean of nsimEIM random variates is taken. Thus nsimEIM controls the accuracy and a larger value may be necessary if the EIMs are not positive-definite. For intercept-only models (y ~ 1) the value of nsimEIM can be smaller (since the common value used is also then taken as the mean over the observations), especially if the number of observations is large.
imethod: An integer with value 1 or 2 or 3 or ... which specifies the initialization method for some parameters or a specific parameter. If failure to converge occurs try the next higher value, and continue until success. For example, imethod = 1 might be the method of moments, and imethod = 2 might be another method. If no value of imethod works then it will be necessary to use arguments such as isigma. For many VGAM family functions it is advisable to try this argument with all possible values to safeguard against problems such as converging to a local solution. VGAM family functions with this argument usually correspond to a model or distribution that is relatively hard to fit successfully, therefore care is needed to ensure the global solution is obtained. So using all possible values that this argument supplies is a good idea.
type.fitted: Character. Type of fitted value returned by the fitted() methods function. The first choice is always the default. The available choices depends on what kind of family function it is. Using the first few letters of the chosen choice is okay. See fittedvlm for more details.
percentiles: Numeric vector, with values between 0 and 100 (although it is not recommended that exactly 0 or 100 be inputted). Used only if type.fitted = "quantiles" or type.fitted = "percentiles", then this argument specifies the values of these quantiles. The argument name tries to reinforce that the values lie between 0 and 100. See fittedvlm for more details.
probs.x, probs.y: Numeric, with values in (0, 1). The probabilites that define quantiles with respect to some vector, usually an x or y of some sort. This is used to create two subsets of data corresponding to `low' and `high' values of x or y. Each value is separately fed into the probs argument of quantile. If the data set size is small then it may be necessary to increase/decrease slightly the first/second values respectively.
lss: Logical. This stands for the ordering: location, scale and shape. Should the ordering of the parameters be in this order? Almost all VGAM family functions have this order by default, but in order to match the arguments of existing R functions, one might need to set lss = FALSE. For example, the arguments of weibullR are scale and shape, whereas rweibull are shape and scale. As a temporary measure (from VGAM 0.9-7 onwards but prior to version 1.0-0), some family functions such as sinmad have an lss argument without a default. For these, setting lss = FALSE will work. Later, lss = TRUE will be the default. Be careful for the dpqr-type functions, e.g., rsinmad.
whitespace: Logical. Should white spaces (" ") be used in the labelling of the linear/additive predictors? Setting TRUE usually results in more readability but it occupies more columns of the output.
oim: Logical. Should the observed information matrices (OIMs) be used for the working weights? In general, setting oim = TRUE means the Newton-Raphson algorithm, and oim = FALSE means Fisher-scoring. The latter uses the EIM, and is usually recommended. If oim = TRUE then nsimEIM is ignored.
nrfs: Numeric, a value in \([0, 1]\). Experimental argument for allowing a mixture of Newton-Raphson and Fisher scoring. The working weights are taken as a linear combination of the two. If nrfs = 0 then Newton-Raphson is used, i.e., the OIM is wholly used. If nrfs = 1 then Fisher scoring is used, i.e., the EIM is wholly used. If convergence is successful then one might expect Newton-Raphson to be faster than Fisher scoring because the former has an order-2 convergence rate while the latter has a linear rate.

zero

Either an integer vector, or a vector of character strings.

If an integer, then it specifies which linear/additive predictor is modelled as intercept-only. That is, the regression coefficients are set to zero for all covariates except for the intercept. If zero is specified then it may be a vector with values from the set \(\{1,2,\ldots,M\}\). The value zero = NULL means model all linear/additive predictors as functions of the explanatory variables. Here, \(M\) is the number of linear/additive predictors. Technically, if zero contains the value \(j\) then the \(j\)th row of every constraint matrix (except for the intercept) consists of all 0 values.

Some VGAM family functions allow the zero argument to accept negative values; if so then its absolute value is recycled over each (usual) response. For example, zero = -2 for the two-parameter negative binomial distribution would mean, for each response, the second linear/additive predictor is modelled as intercepts-only. That is, for all the \(k\) parameters in negbinomial (this VGAM family function can handle a matrix of responses).

Suppose zero = zerovec where zerovec is a vector of negative values. If \(G\) is the usual \(M\) value for a univariate response then the actual values for argument zero are all values in c(abs(zerovec), G + abs(zerovec), 2*G + abs(zerovec), ... ) lying in the integer range \(1\) to \(M\). For example, setting zero = -c(2, 3) for a matrix response of 4 columns with zinegbinomial (which usually has \(G = M = 3\) for a univariate response) would be equivalent to zero = c(2, 3, 5, 6, 8, 9, 11, 12). This example has \(M = 12\). Note that if zerovec contains negative values then their absolute values should be elements from the set 1:G.

Note: zero may have positive and negative values, for example, setting zero = c(-2, 3) in the above example would be equivalent to zero = c(2, 3, 5, 8, 11).

The argument zero also accepts a character vector (for VGAM 1.0-1 onwards). Each value is fed into grep with fixed = TRUE, meaning that wildcards "*" are not useful. See the example below---all the variants work; those with LOCAT issue a warning that that value is unmatched. Importantly, the parameter names are c("location1", "scale1", "location2", "scale2") because there are 2 responses. Yee (2015) described zero for only numerical input. Allowing character input is particularly important when the number of parameters cannot be determined without having the actual data first. For example, with time series data, an ARMA(\(p\),\(q\)) process might have parameters \(\theta_1,\ldots,\theta_p\) which should be intercept-only by default. Then specifying a numerical default value for zero would be too difficult (there are the drift and scale parameters too). However, it is possible with the character representation: zero = "theta" would achieve this. In the future, most VGAM family functions might be converted to the character representation---the advantage being that it is more readable. When programming a VGAM family function that allows character input, the variable predictors.names must be assigned correctly.

If the constraints argument is used then the family function's zero argument (if it exists) needs to be set to NULL. This avoids what could be a probable contradiction. Sometimes setting other arguments related to constraint matrices to FALSE is also a good idea, e.g., parallel = FALSE, exchangeable = FALSE.

ishrinkage

Shrinkage factor \(s\) used for obtaining initial values. Numeric, between 0 and 1. In general, the formula used is something like \(s \mu + (1-s) y\) where \(\mu\) is a measure of central tendency such as a weighted mean or median, and \(y\) is the response vector. For example, the initial values are slight perturbations of the mean towards the actual data. For many types of models this method seems to work well and is often reasonably robust to outliers in the response. Often this argument is only used if the argument imethod is assigned a certain value.

nointercept

An integer-valued vector specifying which linear/additive predictors have no intercepts. Any values must be from the set {1,2,...,\(M\)}. A value of NULL means no such constraints.

multiple.responses

Logical. Some VGAM family functions allow a multivariate or vector response. If so, then usually the response is a matrix with columns corresponding to the individual response variables. They are all fitted simultaneously. Arguments such as parallel may then be useful to allow for relationships between the regressions of each response variable. If multiple.responses = TRUE then sometimes the response is interpreted differently, e.g., posbinomial chooses the first column of a matrix response as success and combines the other columns as failure, but when multiple.responses = TRUE then each column of the response matrix is the number of successes and the weights argument is of the same dimension as the response and contains the number of trials.

earg.link

This argument should be generally ignored.

byrow.arg

Logical. Some VGAM family functions that handle multiple responses have arguments that allow input to be fed in which affect all the responses, e.g., imu for initalizing a mu parameter. In such cases it is sometime more convenient to input one value per response by setting byrow.arg = TRUE; then values are recycled in order to form a matrix of the appropriate dimension. This argument matches byrow in matrix; in fact it is fed into such using matrix(..., byrow = byrow.arg). This argument has no effect when there is one response.

bred

Logical. Some VGAM family functions will allow bias-reduction based on the work by Kosmidis and Firth. Sometimes half-stepping is a good idea; set stepsize = 0.5 and monitor convergence by setting trace = TRUE.

Warning

The zero argument is supplied for convenience but conflicts can arise with other arguments, e.g., the constraints argument of vglm and vgam. See Example 5 below for an example. If not sure, use, e.g., constraints(fit) and coef(fit, matrix = TRUE) to check the result of a fit fit.

The arguments zero and nointercept can be inputted with values that fail. For example, multinomial(zero = 2, nointercept = 1:3) means the second linear/additive predictor is identically zero, which will cause a failure.

Be careful about the use of other potentially contradictory constraints, e.g., multinomial(zero = 2, parallel = TRUE ~ x3). If in doubt, apply constraints() to the fitted object to check.

VGAM family functions with the nsimEIM may have inaccurate working weight matrices. If so, then the standard errors of the regression coefficients may be inaccurate. Thus output from summary(fit), vcov(fit), etc. may be misleading.

Changes relating to the codelss argument have very important consequences and users must beware. Good programming style is to rely on the argument names and not on the order.

Author

T. W. Yee

Details

Full details will be given in documentation yet to be written, at a later date!

References

Yee, T. W. (2015). Vector Generalized Linear and Additive Models: With an Implementation in R. New York, USA: Springer.

Kosmidis, I. and Firth, D. (2009). Bias reduction in exponential family nonlinear models. Biometrika, 96, 793--804.

Miranda-Soberanis, V. F. and Yee, T. W. (2019). New link functions for distribution--specific quantile regression based on vector generalized linear and additive models. Journal of Probability and Statistics, 5, 1--11.

Examples

Run this code

# Example 1
cumulative()
cumulative(link = "probitlink", reverse = TRUE, parallel = TRUE)

# Example 2
wdata <- data.frame(x2 = runif(nn <- 1000))
wdata <- transform(wdata,
         y = rweibull(nn, shape = 2 + exp(1 + x2), scale = exp(-0.5)))
fit <- vglm(y ~ x2, weibullR(lshape = logofflink(offset = -2), zero = 2),
            data = wdata)
coef(fit, mat = TRUE)

# Example 3; multivariate (multiple) response
if (FALSE) {
ndata <- data.frame(x = runif(nn <- 500))
ndata <- transform(ndata,
           y1 = rnbinom(nn, mu = exp(3+x), size = exp(1)),  # k is size
           y2 = rnbinom(nn, mu = exp(2-x), size = exp(0)))
fit <- vglm(cbind(y1, y2) ~ x, negbinomial(zero = -2), data = ndata)
coef(fit, matrix = TRUE)
}
# Example 4
if (FALSE) {
# fit1 and fit2 are equivalent
fit1 <- vglm(ymatrix ~ x2 + x3 + x4 + x5,
             cumulative(parallel = FALSE ~ 1 + x3 + x5), data = cdata)
fit2 <- vglm(ymatrix ~ x2 + x3 + x4 + x5,
             cumulative(parallel = TRUE ~ x2 + x4), data = cdata)
}

# Example 5
udata <- data.frame(x2 = rnorm(nn <- 200))
udata <- transform(udata,
           y1 = rnorm(nn, mean = 1 - 3*x2, sd = exp(1 + 0.2*x2)),
           y2 = rnorm(nn, mean = 1 - 3*x2, sd = exp(1)))
args(uninormal)
fit1 <- vglm(y1 ~ x2, uninormal, data = udata)            # This is okay
fit2 <- vglm(y2 ~ x2, uninormal(zero = 2), data = udata)  # This is okay

# This creates potential conflict
clist <- list("(Intercept)" = diag(2), "x2" = diag(2))
fit3 <- vglm(y2 ~ x2, uninormal(zero = 2), data = udata,
             constraints = clist)  # Conflict!
coef(fit3, matrix = TRUE)  # Shows that clist[["x2"]] was overwritten,
constraints(fit3)  # i.e., 'zero' seems to override the 'constraints' arg

# Example 6 ('whitespace' argument)
pneumo <- transform(pneumo, let = log(exposure.time))
fit1 <- vglm(cbind(normal, mild, severe) ~ let,
             sratio(whitespace = FALSE, parallel = TRUE), data = pneumo)
fit2 <- vglm(cbind(normal, mild, severe) ~ let,
             sratio(whitespace = TRUE,  parallel = TRUE), data = pneumo)
head(predict(fit1), 2)  # No white spaces
head(predict(fit2), 2)  # Uses white spaces

# Example 7 ('zero' argument with character input)
set.seed(123); n <- 1000
ldata <- data.frame(x2 = runif(n))
ldata <- transform(ldata, y1 = rlogis(n, loc = 5*x2, scale = exp(2)))
ldata <- transform(ldata, y2 = rlogis(n, loc = 5*x2, scale = exp(1*x2)))
ldata <- transform(ldata, w1 = runif(n))
ldata <- transform(ldata, w2 = runif(n))
fit7 <- vglm(cbind(y1, y2) ~ x2,
#        logistic(zero = "location1"),  # location1 is intercept-only
#        logistic(zero = "location2"),
#        logistic(zero = "location*"),  # Not okay... all is unmatched
#        logistic(zero = "scale1"),
#        logistic(zero = "scale2"),
#        logistic(zero = "scale"),  # Both scale parameters are matched
         logistic(zero = c("location", "scale2")),  # All but scale1
#        logistic(zero = c("LOCAT", "scale2")),  # Only scale2 is matched
#        logistic(zero = c("LOCAT")),  # Nothing is matched
#        trace = TRUE,
#        weights = cbind(w1, w2),
         weights = w1,
         data = ldata)
coef(fit7, matrix = TRUE)

Run the code above in your browser using DataLab