Description Usage Arguments Details Value Note Author(s) References Examples
stdCoxph
performs regression standardization in Cox proportional hazards models,
at specified values of the exposure, over the sample covariate distribution.
Let T, X, and Z be the survival outcome, the exposure, and a
vector of covariates, respectively. stdCoxph
uses a fitted Cox
proportional hazards model to estimate the standardized
survival function θ(t,x)=E\{S(t|X=x,Z)\}, where t is a specific value of T,
x is a specific value of X, and the expectation is over the marginal
distribution of Z.
1 |
fit |
an object of class |
data |
a data frame containing the variables in the model. This should be the same
data frame as was used to fit the model in |
X |
a string containing the name of the exposure variable X in |
x |
an optional vector containing the specific values of X at which to estimate
the standardized survival function. If X is binary (0/1) or
a factor, then |
t |
an optional vector containing the specific values of T at which to estimate
the standardized survival function. It defaults to all the observed event times
in |
clusterid |
an optional string containing the name of a cluster identification variable when data are clustered. |
subsetnew |
an optional logical statement specifying a subset of observations to be used in the standardization. This set is assumed to be a subset of the subset (if any) that was used to fit the regression model. |
stdCoxph
assumes that a Cox proportional hazards model
λ(t|X,Z)=λ_0(t)exp\{h(X,Z;β)\}
has been fitted. Breslow's estimator of the cumulative baseline hazard Λ_0(t)=\int_0^tλ_0(u)du is used together with the partial likelihood estimate of β to obtain estimates of the survival function S(t|X=x,Z):
\hat{S}(t|X=x,Z)=exp[-\hat{Λ}_0(t)exp\{h(X=x,Z;\hat{β})\}].
For each t in the t
argument and for each x in the x
argument,
these estimates are averaged across all subjects (i.e. all observed values of Z)
to produce estimates
\hat{θ}(t,x)=∑_{i=1}^n \hat{S}(t|X=x,Z_i)/n,
where Z_i is the value of Z for subject i, i=1,...,n. The variance for \hat{θ}(t,x) is obtained by the sandwich formula.
An object of class "stdCoxph"
is a list containing
call |
the matched call. |
input |
|
est |
a matrix with |
vcov |
a list with |
Standardized survival functions are sometimes referred to as (direct) adjusted survival functions in the literature.
stdCoxph
does not currently handle time-varying exposures or covariates.
stdCoxph
internally loops over all values in the t
argument. Therefore,
the function will usually be considerably faster if length(t)
is small.
The variance calculation performed by stdCoxph
does not condition on
the observed covariates \bar{Z}=(Z_1,...,Z_n). To see how this matters,
note that
var\{\hat{θ}(t,x)\}=E[var\{\hat{θ}(t,x)|\bar{Z}\}]+var[E\{\hat{θ}(t,x)|\bar{Z}\}].
The usual parameter β in a Cox proportional hazards model does not depend on \bar{Z}. Thus, E(\hat{β}|\bar{Z}) is independent of \bar{Z} as well (since E(\hat{β}|\bar{Z})=β), so that the term var[E\{\hat{β}|\bar{Z}\}] in the corresponding variance decomposition for var(\hat{β}) becomes equal to 0. However, θ(t,x) depends on \bar{Z} through the average over the sample distribution for Z, and thus the term var[E\{\hat{θ}(t,x)|\bar{Z}\}] is not 0, unless one conditions on \bar{Z}. The variance calculation by Gail and Byar (1986) ignores this term, and thus effectively conditions on \bar{Z}.
Arvid Sjolander
Chang I.M., Gelman G., Pagano M. (1982). Corrected group prognostic curves and summary statistics. Journal of Chronic Diseases 35, 669-674.
Gail M.H. and Byar D.P. (1986). Variance calculations for direct adjusted survival curves, with applications to testing for no treatment effect. Biometrical Journal 28(5), 587-599.
Makuch R.W. (1982). Adjusted survival curve estimation using covariates. Journal of Chronic Diseases 35, 437-443.
Sjolander A. (2016). Regression standardization with the R-package stdReg. European Journal of Epidemiology 31(6), 563-574.
Sjolander A. (2016). Estimation of causal effect measures with the R-package stdReg. European Journal of Epidemiology 33(9), 847-858.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 | require(survival)
n <- 1000
Z <- rnorm(n)
X <- rnorm(n, mean=Z)
T <- rexp(n, rate=exp(X+Z+X*Z)) #survival time
C <- rexp(n, rate=exp(X+Z+X*Z)) #censoring time
U <- pmin(T, C) #time at risk
D <- as.numeric(T < C) #event indicator
dd <- data.frame(Z, X, U, D)
fit <- coxph(formula=Surv(U, D)~X+Z+X*Z, data=dd, method="breslow")
fit.std <- stdCoxph(fit=fit, data=dd, X="X", x=seq(-1,1,0.5), t=1:5)
print(summary(fit.std, t=3))
plot(fit.std)
|
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.