Fit proportional hazards regression models using Q-estimation.
Qcoxph(formula, weights, start, data, knots, wtau = NULL,
control = Qcoxph.control(), ...)an object of classes “Qcoxph”, “coxph”, and “Qest”. See coxph.object for details. All the S3 methods that are available for “coxph” objects will also work with a “Qcoxph” object.
An object of class “Qcoxph” is a list containing at least the following components:
a named vector of coefficients.
the covariance matrix of the coefficients.
number of iterations used.
the vector of linear predictors, one per subject. Note that this vector has not been centered, see predict.coxph for details.
the martingale residuals.
vector of column means of the X matrix. Subsequent survival curves are adjusted to this value.
the number of observations used in the fit.
the number of events used in the fit.
a vector of length 6, containing the number of pairs that are concordant, discordant, tied on x, tied on y, and tied on both, followed by the standard error of the concordance statistic.
other objects used for prediction.
the objective function of the model. Please, interpret with care: read the note in the documentation of Qest.
internal objects.
an object of class “formula” (or one that can be coerced to that class): a symbolic description of the model to be fitted. Use Surv(time, event) ~ x, if the data are right-censored, and Surv(time, time2, event) ~ x, if the data are right-censored and left-truncated (time < time2, time can be -Inf).
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.
optional starting values for the coefficients of the linear predictor.
an optional data frame, list or environment (or object coercible by as.data.frame to a data frame) containing the variables in the model. If not found in data, the variables are taken from environment(formula), typically the environment from which Qcoxph is called.
knots to create the basis of a piecewise linear function. If knots is a vector of at least two elements, it is used to identify the exact position of all knots, including boundaries. If knots is a scalar, its value is used to determine the number of internal knots (knots = 0 is allowed, and fits an Exponential model). If knots is missing, by default max(1, min(floor(n.events/30), 3)) internal knots are used. Unless a vector of knots is provided by the user, the “optimal” position of the knots will be identified using the method described in Muggeo (2008). If this fails, the knots are positioned at the empirical quantiles of the observed events.
an optional function that assigns a different weight to each quantile. By default, all quantiles in (0,1) have the same weight. Please check the documentation of wtrunc for built-in weighting functions.
a list of operational parameters. This is usually passed through Qcoxph.control.
additional arguments for wtau.
Paolo Frumento <paolo.frumento@unipi.it>, Gianluca Sottile <gianluca.sottile@unipa.it>
This function estimates a proportional hazards model, allowing for right-censored and left-truncated data. The syntax and output of Qcoxph are almost identical to those of coxph, but the parameters are estimated using Q-estimation (Sottile and Frumento, 2020). This method can be used to obtain outlier-robust estimators of the regression coefficients.
The quantile function of a proportional hazards model is given by
$$Q(\tau | x) = H0^{-1}(-exp{-x'\beta} log(1 - \tau))$$
where \(H0\) is the baseline cumulative hazard function. In Qcoxph, \(H0\) is parametrized by a piecewise linear function identified by the provided knots. As the number of knots increases, the baseline hazard becomes arbitrarily flexible.
Estimation is carried out by finding the zeroes of a set of integrals equation. The optional argument wtau permits assigning a different weight to each quantile in (0,1). It is possible to choose wtau to be a discontinuous function (e.g., wtau = function(tau){tau < 0.95}). However, this may occasionally result in poorly estimated of the standard errors.
The estimation algorithm is briefly described in the documentation of Qcoxph.control.
Sottile G, and Frumento P (2022). Robust estimation and regression with parametric quantile functions. Computational Statistics and Data Analysis. <doi:10.1016/j.csda.2022.107471>
Muggeo VMR (2008). Segmented: an R package to fit regression models with broken-line relationships. R News 8/1, 20–25.
Qest, for general Q-estimation, and Qlm, for Q-estimation of linear models.
# A proportional-hazards Weibull model
n <- 100
x <- runif(n,0,3)
shape <- 2
t <- rweibull(n, shape = shape, scale = (1/exp(2 + 2*x))^(1/shape)) # time-to-event
c <- runif(n,0,1) # censoring variable
y <- pmin(t,c) # observed response
d <- (t <= c) # event indicator
require(survival)
m1 <- coxph(Surv(y,d) ~ x) # standard Cox model
m2 <- Qcoxph(Surv(y,d) ~ x) # Q-estimator
Run the code above in your browser using DataLab