I want to cross-validate survival curves, and i am not sure I am doing it correctly. The problem is as follows:
I have around 50 patients and their survival times and their (say clinical) measurements. The goal is to use the measurements to find a model that can predict whether a new patient will have low or high risk "to survive". Because of the low number of patients, the goal is just to make sure that the (cross-validated, not single) model can "significantly" distinguish between the low and the high risk group.
For now i did the following: I applied a repeated 10-fold CV, i.e. i repeat 50 times the following procedure:
- Remove the n-th fold from the data.
- Apply the feature selection to the 9 remaining folds [*].
- Apply the training method to the 9 remaining folds to obtain a model.
- Use the model to predict on the n-th fold, i.e. I get a risk for each patient in the n-th fold.
- Gather all the predictions over all folds to have an "independent" prediction on all the training set.
- As the model yields a risk for each patient, I can now take the median of all risks to subdivide the all patients into low and high risk groups. These are now (I hope) cross-validated survival curves. Usually I would apply the log-rank test to get a significance, whether the two curves differ.
But as a single cross-validation is not stable, I did repeat the above procedure 50 times. Then I obtain 50 cross-validated survival curves. My main question is:
- How to compare these 50 curves correctly, using some kind of statistics to yield a "significance"?
One way would be not to generate the cross-validated survival after each repeat, but instead just gather all the 50 (cross-validated) risks for each patient and take the mean of it. With these means i could create two new survival curves and take their log-rank statistics. Which sounds somehow reasonable.
Unluckily i did not find many references for this problem ,the closest is the paper from Simon, Subramanian et al, "Using cross-validation to evaluate predictive accuracy of survival risk classifiers based on high-dimensional data", see https://www.ncbi.nlm.nih.gov/pubmed/21324971.
But it seems that this paper tells me that my approach is not correct, and would be just a "measure" of spread (The paper seems not deal with repeated cross-validation, unluckily).
To obtain a significance (H_0=the groups do not differ at all times) for the two survival groups, the paper tells me I should do a permutation test by repeating these two steps say 1000 times:
- Permute the survival data (time+event) and the variables (this amounts to simply permute the survival data, as my table only consists of (time, event, clinical measurements))
- Apply the full procedure above to get a log-rank statistics.
All the log-rank statistics will now yield a distribution, and i can now determine in which percentile my original, not permuted log-rank statistics will be. Using an alpha level of 5%, this would mean 97.5% percent of all permuted log-rank values should be on the left (or right) of my unpermuted value.
a) Do I have to use the second method (with does not involve repeating cross-validation but permutation) to get a correct result? Or can I work with the repeated cross-validation procedure I alreay have (And, maybe also quite important: Did I correctly understand the procedure in the paper?)
b) Do i have to use the log-rank statistics in the permutation test? Or can I use Harrell's c-index instead?
c) As the 10-fold CV takes several hours (using random forests), is there any kind of shortcut?
d) If my original approach is not correct, can I re-use the 50 repeats i already computed to get a (good) approximation of the "true" p-value?
e) Is there an R package that can help me in any way with the above permutation test?
Thanks a lot.
[*] It would be better to again apply cross-validation in this step, but I ignore this, hoping it will not make a large difference.