VGAM (version 1.1-6)

ordpoisson: Ordinal Poisson Family Function


Fits a Poisson regression where the response is ordinal (the Poisson counts are grouped between known cutpoints).


ordpoisson(cutpoints, countdata = FALSE, NOS = NULL,
           Levels = NULL, = NULL, parallel = FALSE,
           zero = NULL, link = "loglink")



Numeric. The cutpoints, \(K_l\). These must be non-negative integers. Inf values may be included. See below for further details.


Logical. Is the response (LHS of formula) in count-data format? If not then the response is a matrix or vector with values 1, 2, …, L, say, where L is the number of levels. Such input can be generated with cut with argument labels = FALSE. If countdata = TRUE then the response is expected to be in the same format as depvar(fit) where fit is a fitted model with ordpoisson as the VGAM family function. That is, the response is matrix of counts with L columns (if NOS = 1).


Integer. The number of species, or more generally, the number of response random variates. This argument must be specified when countdata = TRUE. Usually NOS = 1.


Integer vector, recycled to length NOS if necessary. The number of levels for each response random variate. This argument should agree with cutpoints. This argument must be specified when countdata = TRUE.

Numeric. Initial values for the means of the Poisson regressions. Recycled to length NOS if necessary. Use this argument if the default initial values fail (the default is to compute an initial value internally).

parallel, zero, link

See poissonff.


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


The input requires care as little to no checking is done. If fit is the fitted object, have a look at fit@extra and depvar(fit) to check.


This VGAM family function uses maximum likelihood estimation (Fisher scoring) to fit a Poisson regression to each column of a matrix response. The data, however, is ordinal, and is obtained from known integer cutpoints. Here, \(l=1,\ldots,L\) where \(L\) (\(L \geq 2\)) is the number of levels. In more detail, let \(Y^*=l\) if \(K_{l-1} < Y \leq K_{l}\) where the \(K_l\) are the cutpoints. We have \(K_0=-\infty\) and \(K_L=\infty\). The response for this family function corresponds to \(Y^*\) but we are really interested in the Poisson regression of \(Y\).

If NOS=1 then the argument cutpoints is a vector \((K_1,K_2,\ldots,K_L)\) where the last value (Inf) is optional. If NOS>1 then the vector should have NOS-1 Inf values separating the cutpoints. For example, if there are NOS=3 responses, then something like ordpoisson(cut = c(0, 5, 10, Inf, 20, 30, Inf, 0, 10, 40, Inf)) is valid.


Yee, T. W. (2020). Ordinal ordination with normalizing link functions for count data, (in preparation).

See Also

poissonff, polf, ordered.


Run this code
set.seed(123)  # Example 1
x2 <- runif(n <- 1000); x3 <- runif(n)
mymu <- exp(3 - 1 * x2 + 2 * x3)
y1 <- rpois(n, lambda = mymu)
cutpts <- c(-Inf, 20, 30, Inf)
fcutpts <- cutpts[is.finite(cutpts)]  # finite cutpoints
ystar <- cut(y1, breaks = cutpts, labels = FALSE)
# }
plot(x2, x3, col = ystar, pch = as.character(ystar))
# }
table(ystar) / sum(table(ystar))
fit <- vglm(ystar ~ x2 + x3, fam = ordpoisson(cutpoi = fcutpts))
head(depvar(fit))  # This can be input if countdata = TRUE
coef(fit, matrix = TRUE)

# Example 2: multivariate and there are no obsns between some cutpoints
cutpts2 <- c(-Inf, 0, 9, 10, 20, 70, 200, 201, Inf)
fcutpts2 <- cutpts2[is.finite(cutpts2)]  # finite cutpoints
y2 <- rpois(n, lambda = mymu)   # Same model as y1
ystar2 <- cut(y2, breaks = cutpts2, labels = FALSE)
table(ystar2) / sum(table(ystar2))
fit <- vglm(cbind(ystar,ystar2) ~ x2 + x3, fam =
            ordpoisson(cutpoi = c(fcutpts,Inf,fcutpts2,Inf),
                       Levels = c(length(fcutpts)+1,length(fcutpts2)+1),
                       parallel = TRUE), trace = TRUE)
coef(fit, matrix = TRUE)
summary(depvar(fit))  # Some columns have all zeros
# }

Run the code above in your browser using DataCamp Workspace