Overview
A basic component for our modelling of multivariate survival data isthat many models are build around marginals that on Cox form. Themarginal Cox model can be fitted efficiently in the mets package, inparticular the handling of strata and robust standard errors isoptimized.
The basic models assumes that each subject has a marginal on Cox-form\[\lambda_{g(k,i)}(t) \exp( X_{ki}^T \beta).\] where\(g(k,i)\) gives thestrata for the subject.
We here discuss and show how to get robust standard errors of
the regression parameters
the baseline
and how to do goodness of fit test using
- cumulative residuals score test
First we generate some data from the Clayton-Oakes model, with\(5\) members in each cluster and a varianceparameter at\(2\)
library(mets)options(warn=-1)set.seed(1000)# to control output in simulatins for p-values below. n<-1000 k<-5 theta<-2 data<-simClaytonOakes(n,k,theta,0.3,3)The data is on has one subject per row.
- time : time of event
- status : 1 for event and 0 for censoring
- x : x is a binary covariate
- cluster : cluster
Now we fit the model and produce robust standard errors for bothregression parameters and baseline.
First, recall that the baseline for strata\(g\) is asymptotically equivalent to\[\begin{align}\hat A_g(t) - A_g(t) & = \sum_{k \in g} \left( \sum_{i: ki \in g}\int_0^t \frac{1}{S_{0,g}} dM_{ki}^g - P^g(t) \beta_k \right)\end{align}\] with\(P^g(t) = \int_0^tE_g(s) d \hat \Lambda_g(s)\) the derivative of\(\int_0^t 1/S_{0,g}(s) dN_{\cdot g}\) wrt to\(\beta\), and\[\begin{align}\hat \beta - \beta & = I(\tau)^{-1} \sum_k ( \sum_i \int_0^\tau(Z_{ki} - E_{g}) dM_{ki}^g ) = \sum_k \beta_{k}\end{align}\] with\[\begin{align}M_{ki}^g(t) & = N_{ki}(t) - \int_0^t Y_{ki}(s) \exp( Z_{ki} \beta)d \Lambda_{g(k,i)}(t), \\\beta_{k} & = I(\tau)^{-1} \sum_i \int_0^\tau (Z_{ki} - E_{g})dM_{ki}^g\end{align}\] the basic 0-mean processes, that are martingales inthe iid setting, and\(I(t)\) is thederivative of the total score,\(\hatU(t,\beta))\), with respect to\(\beta\) evaluated at time\(t\).
The variance of the baseline of strata g is estimated by\[\begin{align}\sum_{k \in g} ( \sum_{i: ki \in g} \int_0^t \frac{1}{S_{0,g(k,i)}}d\hat M_{ki}^g - P^g(t) \beta_k )^2\end{align}\] that can be computed using the particular structureof\[\begin{align}d \hat M_{ik}^g(t) & = dN_{ik}(t) - \frac{1}{S_{0,g(i,k)}}\exp(Z_{ik} \beta) dN_{g.}(t)\end{align}\]
This robust variance of the baseline and the iid decomposition for\(\beta\) is computed in mets as:
out<-phreg(Surv(time,status)~x+cluster(cluster),data=data)summary(out)#>#> n events#> 5000 4854#> coeffients:#> Estimate S.E. dU^-1/2 P-value#> x 0.287859 0.028177 0.028897 0#>#> exp(coeffients):#> Estimate 2.5% 97.5%#> x 1.3336 1.2619 1.4093# robust standard errors attached to output rob<-robust.phreg(out)We can get the iid decomposition of the\(\hat \beta - \beta\) by
# making iid decomposition of regression parameters betaiid<-IC(out)head(betaiid)#> x#> 1 -0.34616008#> 2 -1.44918926#> 3 -0.03898156#> 4 0.42156050#> 5 0.34253904#> 6 -0.07706668# robust standard errorscrossprod(betaiid/NROW(betaiid))^.5#> x#> x 0.02817714# same asWe now look at the plot with robust standard errors
We can also make survival prediction with robust standard errorsusing the phreg.
Finally, just to check that we can recover the model we also estimatethe dependence parameter
tt<-twostageMLE(out,data=data)summary(tt)#> Dependence parameter for Clayton-Oakes model#> Variance of Gamma distributed random effects#> $estimates#> Coef. SE z P-val Kendall tau SE#> dependence1 0.5316753 0.03497789 15.20032 0 0.2100093 0.0109146#>#> $type#> NULL#>#> attr(,"class")#> [1] "summary.mets.twostage"