Fits Ghosh-Lin IPCW Cox-type model
recreg(
formula,
data,
cause = 1,
death.code = 2,
cens.code = 0,
cens.model = ~1,
weights = NULL,
offset = NULL,
Gc = NULL,
wcomp = NULL,
marks = NULL,
augmentation.type = c("lindyn.augment", "lin.augment"),
...
)
formula with 'Event' outcome
data frame
of interest (1 default)
codes for death (terminating event, 2 default)
code of censoring (0 default)
for stratified Cox model without covariates
weights for score equations
offsets for model
censoring weights for time argument, default is to calculate these with a Kaplan-Meier estimator, should then give G_c(T_i-)
weights for composite outcome, so when cause=c(1,3), we might have wcomp=c(1,2).
a mark value can be specified, this is vector from the data-frame where the mark value can be found at all events
of augmentation when augmentation model is given
Additional arguments to lower level funtions
For Cox type model : $$ E(dN_1(t)|X) = \mu_0(t)dt exp(X^T \beta) $$ by solving Cox-type IPCW weighted score equations $$ \int (Z - E(t)) w(t) dN_1(t) $$ where $$w(t) = G(t) (I(T_i \wedge t < C_i)/G_c(T_i \wedge t))$$ and $$E(t) = S_1(t)/S_0(t)$$ and $$S_j(t) = \sum X_i^j w_i(t) \exp(X_i^T \beta)$$.
The iid decomposition of the beta's are on the form $$ \int (Z - E ) w(t) dM_1 + \int q(s)/p(s) dM_c $$ and returned as iid.
Events, deaths and censorings are specified via stop start structure and the Event call, that via a status vector and cause (code), censoring-codes (cens.code) and death-codes (death.code) indentifies these. See example and vignette.
## data with no ties
library(mets)
data(hfactioncpx12)
hf <- hfactioncpx12
hf$x <- as.numeric(hf$treatment)
dd <- data.frame(treatment=levels(hf$treatment),id=1)
gl <- recreg(Event(entry,time,status)~treatment+cluster(id),data=hf,cause=1,death.code=2)
summary(gl)
#>
#> n events
#> 2132 1391
#>
#> 741 clusters
#> coeffients:
#> Estimate S.E. dU^-1/2 P-value
#> treatment1 -0.110499 0.078650 0.053776 0.16
#>
#> exp(coeffients):
#> Estimate 2.5% 97.5%
#> treatment1 0.89539 0.76747 1.0446
#>
#>
pgl <- predict(gl,dd,se=1); plot(pgl,se=1)
## censoring stratified after treatment
gls <- recreg(Event(entry,time,status)~treatment+cluster(id),data=hf,
cause=1,death.code=2,cens.model=~strata(treatment))
summary(gls)
#>
#> n events
#> 2132 1391
#>
#> 741 clusters
#> coeffients:
#> Estimate S.E. dU^-1/2 P-value
#> treatment1 -0.109604 0.078701 0.053777 0.1637
#>
#> exp(coeffients):
#> Estimate 2.5% 97.5%
#> treatment1 0.89619 0.76809 1.0457
#>
#>
## IPCW at 2 years
ll2 <- recregIPCW(Event(entry,time,status)~treatment+cluster(id),data=hf,
cause=1,death.code=2,time=2,cens.model=~strata(treatment))
summary(ll2)
#>
#> n events
#> 741 0
#>
#> 741 clusters
#> coeffients:
#> Estimate Std.Err 2.5% 97.5% P-value
#> (Intercept) 0.452257 0.062479 0.329801 0.574714 0.0000
#> treatment1 -0.078348 0.097213 -0.268883 0.112187 0.4203
#>
#> exp(coeffients):
#> Estimate 2.5% 97.5%
#> (Intercept) 1.57186 1.39069 1.7766
#> treatment1 0.92464 0.76423 1.1187
#>
#>