ctqr (version 2.1)

ctqr: Censored and Truncated Quantile Regression

Description

Fits a quantile regression model to possibly censored and truncated data, e.g., survival data.

Usage

ctqr(formula, data, weights, p = 0.5, CDF, control = ctqr.control(), ...)

Value

An object of class “ctqr”, which is a list with the following items:

p

the quantile(s) being estimated.

coefficients

a named vector or matrix of quantile regression coefficients.

call

the matched call.

n.it

the number of iterations.

converged

logical. The convergence status.

fitted

the fitted values.

terms

the terms object used.

mf

the model frame used.

covar

the estimated asymptotic covariance matrix.

CDF

the used CDF.

Note that the dimension of all items, except call, terms, mf, and CDF, is the same as the dimension of p. For example, if p = c(0.25,0.5,0.75), coefficients

and fitted will be 3-columns matrices; n.it and converged will be vectors of 3 elements; and covar will be a list of three covariance matrices.

The generic accessor functions summary, plot, predict, coef, terms, nobs, can be used to extract information from the model. The functions waldtest (from the package lmtest), and linearHypothesis (from the package car) can be used to perform Wald test, and to test linear restrictions. These functions, however, will only work if p is scalar.

Arguments

formula

an object of class “formula”: a symbolic description of the regression model. The response must be a Surv object as returned by Surv (see ‘Details’).

data

an optional data frame containing the variables in the model.

weights

an optional vector of weights to be used in the fitting process. The weights will always be normalized to sum to the sample size. This implies that, for example, using double weights will not halve the standard errors.

p

numerical vector indicating the order of the quantile(s) to be fitted.

CDF

an object of class “pch”, i.e., the result of a call to pchreg. If missing, it will be computed internally with default settings. See ‘Details’.

control

a list of operational parameters for the optimization algorithm, usually passed via ctqr.control.

...

for future arguments.

Author

Paolo Frumento <paolo.frumento@unipi.it>

Details

This function implements the method described in Frumento and Bottai (2017) for censored, truncated quantile regression, and Frumento (2022) for interval-censored quantile regression.

The left side of formula must be of the form Surv(time, event) if the data are right-censored, Surv(time0, time, event) if the data are right-censored and left-truncated (time0 < time, time0 can be -Inf), and Surv(time1, time2, type = "interval2") if the data are interval-censored (use time1 = time2 for exact observations, time1 = -Inf or NA for left-censored, and time2 = Inf or NA for right-censored). Using Surv(time) is also allowed and indicates that the data are neither censored nor truncated.

The conditional distribution function (CDF) of the response variable represents a nuisance parameter and is estimated preliminarly via pchreg. If missing, CDF = pchreg(formula) is used as default. See the “Note” and the documentation of pchreg.

Estimation is carried out using an algorithm for gradient-based optimization. To estimate the asymptotic covariance matrix, standard two-step procedures are used (e.g., Ackerberg et al., 2012).

References

Ackerberg, D., Chen, X., and Hahn, J. (2012). A practical asymptotic variance estimator for two-step semiparametric estimators. The Review of Economics and Statistics, 94 (2), 481-498.

Frumento, P., and Bottai, M. (2017). An estimating equation for censored and truncated quantile regression. Computational Statistics and Data Analysis, Vol.113, pp.53-63. ISSN: 0167-9473.

Frumento, P. (2022). A quantile regression estimator for interval-censored data. The International Journal of Biostatistics, 19 (1), pp. 81-96.

See Also

plot.ctqr, predict.ctqr, pchreg

Examples

Run this code
# Using simulated data


# Example 1 - censored data ####################################################

n <- 1000
x1 <- runif(n); x2 <- runif(n)   # covariates
t <- runif(n, 0, 1 + x1 + x2)    # time variable (e.g., time to death)
c <- runif(n,0,5)                # censoring variable (e.g., end of follow-up)
y <- pmin(t,c)                   # observed variable = min(t,c)
d <- (t <= c)                    # 1 = event (e.g., death), 0 = censored

CDF1 <- pchreg(Surv(y,d) ~ x1 + x2)
model1 <- ctqr(Surv(y,d) ~ x1 + x2, p = 0.5, CDF = CDF1)
model2 <- ctqr(Surv(y,d) ~ x1, p = 0.5, CDF = CDF1)

# model1 is identical to ctqr(Surv(y,d) ~ x1 + x2, p = 0.5)
# model2 is NOT identical to ctqr(Surv(y,d) ~ x1, p = 0.5), 
  # which would have default CDF = pchreg(Surv(y,d) ~ x1)


# Example 2 - censored and truncated data ######################################

n <- 1000
x1 <- runif(n); x2 <- runif(n)   # covariates
t <- runif(n, 0, 1 + x1 + x2)    # time variable
c <- runif(n,0,5)                # censoring variable
y <- pmin(t,c)                   # observed variable = min(t,c)
d <- (t <= c)                    # 1 = event, 0 = censored

z <- rnorm(n) # truncation variable (e.g., time at enrollment)
w <- which(y > z) # data are only observed when y > z
z <- z[w]; y <- y[w]; d <- d[w]; x1 <- x1[w]; x2 <- x2[w]

# implement various CDFs and choose the model with smallest AIC

CDFs <- list(
  pchreg(Surv(z,y,d) ~ x1 + x2, breaks = 5),
  pchreg(Surv(z,y,d) ~ x1 + x2, breaks = 10),
  pchreg(Surv(z,y,d) ~ x1 + x2 + x1:x2, breaks = 5),
  pchreg(Surv(z,y,d) ~ x1 + x2 + x1^2 + x2^2, breaks = 10)
)

CDF <- CDFs[[which.min(sapply(CDFs, function(obj) AIC(obj)))]]
summary(ctqr(Surv(z,y,d) ~ x1 + x2, p = 0.5, CDF = CDF))


# Example 3 - interval-censored data ###########################################
# t is only known to be in the interval (t1,t2) ################################

n <- 1000
x1 <- runif(n); x2 <- runif(n)     # covariates
t <- runif(n, 0, 10*(1 + x1 + x2)) # time variable
t1 <- floor(t)                     # lower extreme of the interval
t2 <- ceiling(t)                   # upper extreme of the interval

model <- ctqr(Surv(t1,t2, type = "interval2") ~ x1 + x2, p = 0.5)

Run the code above in your browser using DataLab