cenpoisson {VGAM}R Documentation

Censored Poisson Family Function

Description

Family function for a censored Poisson response.

Usage

cenpoisson(link = "loge", earg = list(), imu = NULL)

Arguments

link Link function applied to the mean. See Links for more choices.
earg Extra argument optionally used by the link function. See Links for more information.
imu Optional initial value. See CommonVGAMffArguments for more information.

Details

Often a table of Poisson counts has an entry J+ meaning >= J. This family function is similar to poissonff but handles such censored data. The input requires Surv. Only a univariate response is allowed. The Newton-Raphson algorithm is used.

Value

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

Warning

As the response is discrete, care is required with Surv, especially with "interval" censored data because of the (start, end] format. See the examples below. The examples have y < L as left censored and y >= U (formatted as U+) as right censored observations, therefore L <= y < U is for uncensored and/or interval censored observations. Consequently the input must be tweaked to conform to the (start, end] format.

Note

The function poissonff should be used when there are no censored observations. Also, NAs are not permitted with Surv, nor is type="counting".

Author(s)

Thomas W. Yee

References

See survival for background.

See Also

Surv, poissonff, Links.

Examples

# Example 1: right censored data
set.seed(123)
y = rpois(n <- 100, exp(3))
U = 20
cy = pmin(U, y)
rcensored = y >= U
table(cy)
table(rcensored)
status = ifelse(rcensored, 0, 1)
table(i <- print(Surv(cy, status)))  # Check; U+ means >= U
fit = vglm(Surv(cy, status) ~ 1, cenpoisson, trace=TRUE)
coef(fit, mat=TRUE)
table(print(fit@y))  # Another check; U+ means >= U

# Example 2: left censored data
L = 15
cy = pmax(L, y)
lcensored = y <  L   # Note y < L, not cy == L or y <= L
table(cy)
table(lcensored)
status = ifelse(lcensored, 0, 1)
table(i <- print(Surv(cy, status, type="left")))  # Check
fit = vglm(Surv(cy, status, type="left") ~ 1, cenpoisson, trace=TRUE)
coef(fit, mat=TRUE)

# Example 3: interval censored data
Lvec = rep(L, len=n)
Uvec = rep(U, len=n)
icensored = Lvec <= y & y < Uvec  # Neither lcensored or rcensored
table(icensored)
status = rep(3, n)                    # 3 means interval censored
status = ifelse(rcensored, 0, status) # 0 means right censored
status = ifelse(lcensored, 2, status) # 2 means left  censored
# Have to adjust Lvec and Uvec because of the (start, end] format:
Lvec[icensored] = Lvec[icensored] - 1
Uvec[icensored] = Uvec[icensored] - 1
Lvec[lcensored] = Lvec[lcensored]  # Remains unchanged
Lvec[rcensored] = Uvec[rcensored]  # Remains unchanged
table(i <- print(Surv(Lvec, Uvec, status, type="interval")))  # Check

fit = vglm(Surv(Lvec, Uvec, status, type="interval") ~ 1,
           cenpoisson, trace=TRUE)
coef(fit, mat=TRUE)
table(print(fit@y))  # Another check

# Example 4: Add in some uncensored observations
index = (1:n)[icensored]
index = index[1:min(4,length(index))]
status[index] = 1 # actual or uncensored value
Lvec[index] = y[index]
table(i <- print(Surv(Lvec, Uvec, status, type="interval")))  # Check

fit = vglm(Surv(Lvec, Uvec, status, type="interval") ~ 1,
           cenpoisson, trace=TRUE, crit="c")
coef(fit, mat=TRUE)
table(print(fit@y))  # Another check

[Package VGAM version 0.7-7 Index]