Is it possible to create survfit and survdiff objects using imputed data?

Question

I have a dataset in which 45% of participants have missing data. Given the high proportion of missingness, I planned to conduct survival analyses on both an imputed and non-imputed dataset.

I am confident the data are MAR and have used the mice package in R to impute the missing values. I have found it easy enough to create a pooled Cox regression model, using the with() and pool() functions.

I would now like to create survfit and survdiff objects for the construction of Kaplan-Meier curves, but cannot find a way to perform this in R. As an example, if I use the following code:

pathsurv <- with(data = imputed_adenosurv,survfit (Surv (newos, newcensored) ~ pathology))

pathsurv <- with (data = imputed_adenosurv, survdiff(Surv(newos, newcensored) ~ pathology))

I simply get a list with the results of 45 different imputations. Using the pool() function does not rectify this.

Would appreciate any help that can be offered :)

Answer 1

The survdiff() and survfit() functions don't seem to return estimates of coefficients and coefficient-covariance matrices in a form that the pool() function in mice can readily use. You can apply Rubin's rules to get pooled estimates including variances, following the principles in Stef van Buuren's book.

survdiff()

An object returned by the survdiff() function with the default log-rank test contains the numbers of observed and expected (under the null hypothesis) events for each stratum, in the obs and exp slots respectively. The var slot contains the corresponding covariance matrix. If there are N strata there are only N-1 independent contributions to obs-exp, so one stratum needs to be removed for the statistical test, which is a multi-parameter Wald test of the null hypothesis that all of the remaining obs-exp values are 0.

Section 5.3.1 of van Buuren's book shows how to pool such estimates. A Wald statistic is calculated from the imputation-adjusted pooled covariance matrix (see Section 2.3.2), then evaluated in an F-test having degrees of freedom based on the number of imputations, the number of coefficients, and the fraction of missing information.

The D1() function in mice implements this. I'm not sure how to get the survdiff() output into the needed form. "Any R expression produced by expression() can be evaluated on the multiply imputed data," however (Section 5.1.4). One might write an expression to put the needed (linearly independent) obs-exp and covariance estimates from the survdiff() models into a format that the broom package can interpret to extract model components for pool(). (See the documentation for pool().) I don't have any experience with that, however.

survfit()

The pooled Kaplan-Meier curves would be the averages of the individual curves, by stratum. The object returned by survfit() contains standard errors for each event time point, of the type that was specified in the function call. To get closer to the normality assumed by Rubin's rules, choose the log (default) or the log-log type, apply Rubin's rules to pool at each time point, then re-transform back to the survival scale.