plot.survival.rfsrc: Survival Prediction Performance and Diagnostic Plots

In randomForestSRC: Fast Unified Random Forests for Survival, Regression, and Classification (RF-SRC)

Survival Prediction Performance and Diagnostic Plots

Description

Plot survival estimates and calculate inverse-probability-of-censoring weighted Brier scores and cumulative/dynamic time-dependent AUC. plotBrierAUC also compares the forest Brier curve with a predictor-free Kaplan–Meier reference estimated from the training data. get.cindex calculates concordance-based prediction error from observed outcomes and risk predictions for right-censored survival or competing risks.

Usage

## S3 method for class 'rfsrc'
plot.survival(x, show.plots = TRUE, subset,
  collapse = FALSE, cens.model = c("km", "rfsrc"), ...)

get.brier.survival(o, subset = NULL,
  cens.model = c("km", "rfsrc"), papply = lapply,
  times = NULL, conf.int = FALSE, keep.matrix = TRUE)

get.auct.survival(o, subset = NULL,
  cens.model = c("km", "rfsrc"), papply = lapply,
  times = NULL, conf.int = FALSE)

plotBrierAUC(x, subset = NULL,
  cens.model = c("km", "rfsrc"), papply = lapply,
  times = NULL, conf.int = TRUE,
  plots = c("brier", "auct"), show.plots = TRUE,
  brier.null = TRUE, ...)

get.cindex(time, censoring, predicted, weight, fast, do.trace = FALSE)

Arguments

x, o	An object of class (rfsrc, grow) or (rfsrc, predict). The two numerical helpers and plotBrierAUC also accept an object of class (rfsrc, forest).
show.plots	Should plots be displayed?
subset	Vector indicating which predicted cases are to be used. All cases are used if not specified.
collapse	Collapse the individual survival functions into their ensemble mean?
cens.model	Method used to estimate the censoring distribution for inverse probability of censoring weighting. In all cases the censoring model is estimated from the full grow data, not from the evaluation outcomes: km: Kaplan–Meier estimator. rfsrc: Random survival forest estimator of the conditional censoring distribution. When censoring is present, this option requires the grow and evaluation covariates to be stored in the object. Because rfsrc.fast uses forest=FALSE by default, use cens.model="km" for its default reduced object or refit with forest=TRUE.
papply	Function used in place of lapply for repeated calculations.
times	Optional numeric vector of prediction horizons. The forest survival curves are evaluated as right-continuous step functions. The forest time grid is used by default.
conf.int	Controls pointwise confidence intervals. Use FALSE for no intervals, TRUE for 95 percent intervals, or a numeric value strictly between zero and one for a different confidence level.
keep.matrix	Should the subject-by-time matrix of IPCW Brier contributions be retained? This matrix is used by plot.survival for mortality-stratified curves.
plots	One or both of "brier" and "auct". "auc" is also accepted as an alias for "auct".
brier.null	For plotBrierAUC, calculate and display the Brier score of a predictor-free Kaplan–Meier survival estimate from the full grow data? The default is TRUE. The reference is also returned when show.plots=FALSE. This argument has no effect when only AUC is requested.
...	Further graphical arguments. For plotBrierAUC, xlim, ylim, titles, and axis options apply to each requested panel. The arguments col, lwd, lty, and type control the forest curve; band.col and band.border control its confidence band. Use null.col, null.lwd, and null.lty to style the Brier reference (defaults: black, width 2, dashed). Use legend.pos to position its Forest/Null legend (default "topright"), or legend.pos=FALSE to omit the legend.
time	For get.cindex, the observed follow-up times, aligned with censoring and the prediction rows.
censoring	For get.cindex, numeric status codes: zero for censoring and one for an event in right-censored survival. For competing risks, use event codes 1, ..., J.
predicted	For get.cindex, risk predictions, with larger values representing higher risk. For right-censored survival, a vector or one-column matrix or data frame. For competing risks, a matrix or data frame with column j corresponding to event code j, or a list of event-specific prediction vectors. The first J columns or list entries are used.
weight	Optional numeric vector of precomputed concordance weights, one per observation. Omit it for the unweighted calculation. The helper does not estimate censoring weights from the supplied outcomes.
fast	Optional native concordance-algorithm selector. Normally leave it unspecified so that the native library chooses the implementation. It is not used by the weighted competing-risk branch.
do.trace	For get.cindex, enable native diagnostic output?

Details

plot.survival produces the following plots, going from top to bottom and left to right:

1. Forest estimated survival function for each individual. The thick red line is the ensemble mean survival and the thick green line is the marginal Nelson–Aalen survival estimate.

2. Brier score stratified by ensemble mortality. Stratification is into four groups corresponding to the 0–25, 25–50, 50–75 and 75–100 percentile ranges of mortality. The red line is the overall Brier score.

3. Cumulative integrated Brier score divided by elapsed time, labeled CRPS in the plot.

4. Mortality (Ishwaran et al., 2008) versus observed time. Blue points are events and black points are censored observations. Mortality is estimated risk calibrated to the scale of the number of events. For example, a mortality value of 100 means that if all individuals had the same covariate values, an average of 100 events would be expected.

Whenever possible, out-of-bag predictions are used. For a prediction object, the predicted test survival curves and test outcomes are used for scoring, while the censoring distribution is still estimated from the grow data. If a prediction object contains no outcomes, plot.survival displays only its predicted survival curves. Brier score, AUC, CRPS, and mortality-versus-time cannot be calculated without evaluation outcomes. Direct calls to get.brier.survival, get.auct.survival, or plotBrierAUC therefore require outcomes for the cases being evaluated.

For subject i at horizon t, the IPCW Brier contribution is

L_i(t) = \frac{I(T_i \le t, \Delta_i > 0)}{\widehat G(T_i-)} \widehat S_i(t)^2 + \frac{I(T_i > t)}{\widehat G(t)} \{1-\widehat S_i(t)\}^2,

with subject-specific values of \widehat G when cens.model="rfsrc". The Brier score is the sample mean of these contributions.

The time-dependent AUC uses a cumulative/dynamic definition. Cases at t are observed events satisfying T_i \le t; controls satisfy T_i > t; observations censored at or before t are not included in the case-control comparison. The time-specific risk score is 1-\widehat S_i(t), and comparisons are weighted by the same grow-data censoring distribution.

When confidence intervals are requested, the Brier standard error is obtained by deleting one subject-level IPCW loss at a time while holding the survival predictions and censoring weights fixed.

The AUC standard error uses a stratified delete-one jackknife. One observed case or control is deleted at a time, the remaining IPCW weights within that stratum are renormalized, and the resulting case-delete and control-delete variance components are added. With equal weights this calculation reduces exactly to the ordinary DeLong variance. The effective case and control sample sizes and the largest normalized IPCW weights are returned to diagnose horizons at which a few observations dominate the weighted comparison.

Both are pointwise conditional standard errors, and the reported intervals are pointwise normal approximations rather than simultaneous bands. The survival predictions and the estimated censoring distribution are treated as fixed. They do not account for uncertainty from training the survival forest or estimating the censoring model. For an independent prediction sample, the intervals have a direct conditional test-performance interpretation. For a grow object, shared out-of-bag fits and the shared censoring estimate induce dependence among subject-level quantities, so the intervals should be viewed as fixed-fit working intervals rather than repeated-training confidence intervals. A score and its interval are returned as NA at a horizon where a required censoring survival probability is zero. The AUC standard error is also NA if fewer than two cases or two controls are available, or if deleting one subject leaves no positive IPCW weight in its stratum.

plotBrierAUC is a base-R plotting helper. It calls a shared numerical engine so that Brier score, its null reference, and AUC use the same censoring model, then adds the requested pointwise confidence bands for the forest curves. The AUC panel retains its horizontal reference at 0.5.

Null reference for the Brier score

With brier.null=TRUE, the Brier panel includes a dashed predictor-free reference and a Forest/Null legend. Its survival prediction is the same for every evaluated subject:

\widehat S_0(t) = \prod_{u \le t}\left\{1-\frac{d_{\mathrm{grow}}(u)} {Y_{\mathrm{grow}}(u)}\right\},

where d_{\mathrm{grow}}(u) and Y_{\mathrm{grow}}(u) are the event count and risk-set size at grow-data event time u. All finite grow outcomes are used, including their frequencies at tied times. Events are processed before censorings at a tied time. This is the event-survival estimate; cens.model separately controls the censoring distribution used for scoring.

The null Brier curve is obtained by replacing \widehat S_i(t) with \widehat S_0(t) in the IPCW loss above. It uses the same horizons, evaluation subset, censoring weights, and finite forest-loss contributions as the forest curve. Censored observations at or before a horizon retain their zero contribution in both averages. A reference horizon is NA if the forest horizon is unsupported or a required null loss is nonfinite. Confidence bands are drawn only for the forest.

For test evaluation, the reference is estimated from the saved grow outcomes, not from the test outcomes. Changing subset changes the cases being scored, not the grow-data reference. For training evaluation, the reference uses the full training sample and is an in-sample benchmark, even when the forest predictions are OOB. The last Kaplan–Meier estimate is carried forward after its last observed time, subject to the existing censoring-support checks.

A forest curve below the null curve has lower estimated prediction error than this predictor-free reference at that horizon. The two curves are compared directly on the Brier scale; the score is not rescaled and no comparison test is added.

Only right-censored survival families are supported. Competing-risk analyses should use plot.competing.risk.

Concordance error

get.cindex accepts aligned outcome and risk-prediction values directly. Supply predicted.oob for OOB evaluation or predicted with the corresponding outcomes for new-data evaluation. It uses the supplied values without fitting a forest or selecting an OOB component automatically.

A status code greater than one selects the competing-risk branch. For that branch, J is the largest usable status code and prediction columns are indexed by event code, not matched by their names. Without weights, the calculation for event j uses observations censored or experiencing event j; observations with other events are excluded. With weights, it uses the native event-specific weighted concordance calculation instead.

Missing times, statuses, predictions, and supplied weights are excluded from the applicable calculation. The helper returns the native concordance error; it is distinct from the time-specific cumulative/dynamic AUC returned by get.auct.survival.

Value

get.brier.survival returns a list containing brier.score, the optional brier.matx, integrated scores crps and crps.std, the estimated censoring distribution, aligned grow and evaluation event information, survival predictions, mortality, and the selected subset. When confidence intervals are requested, brier.score also contains std.err, lower, upper, and n.eval.

get.auct.survival returns a corresponding list whose auct.score contains time, auct, n.case, n.control, effective sample sizes n.case.eff and n.control.eff, and largest normalized weights max.case.weight and max.control.weight. The columns std.err, lower, and upper are included when confidence intervals are requested.

plotBrierAUC invisibly returns a list containing the requested Brier and/or AUC result objects. When Brier is requested with brier.null=TRUE, the Brier result also contains a null list with brier.score (columns time and brier.score), integrated scores crps and crps.std, and the number of evaluated contributions n.eval. It also contains the common null survival vector on the scoring grid, n.train (finite grow outcomes), method="kaplan-meier", and source="grow". No subject-by-time null loss matrix or null confidence band is returned. The integrated null scores use the same trapezoidal rule and standardization as the forest scores. Direct numerical-helper calls retain their existing return structure.

With evaluation outcomes, plot.survival invisibly returns the mortality-stratified and overall Brier and cumulative integrated Brier curves. Without evaluation outcomes, it invisibly returns the predicted survival curves and their ensemble mean.

get.cindex returns a numeric concordance-error value for right-censored survival, or a numeric vector in event-code order 1, ..., J for competing risks. In the unweighted competing-risk calculation an event with fewer than two usable observations returns NA.

Author(s)

Hemant Ishwaran and Udaya B. Kogalur

References

Gerds T.A. and Schumacher M. (2006). Consistent estimation of the expected Brier score in general survival models with right-censored event times, Biometrical Journal, 48:1029–1040.

Graf E., Schmoor C., Sauerbrei W. and Schumacher M. (1999). Assessment and comparison of prognostic classification schemes for survival data, Statistics in Medicine, 18:2529–2545.

DeLong E.R., DeLong D.M. and Clarke-Pearson D.L. (1988). Comparing the areas under two or more correlated receiver operating characteristic curves: a nonparametric approach, Biometrics, 44:837–845.

Efron B. and Tibshirani R.J. (1993). An Introduction to the Bootstrap. Chapman and Hall, New York.

Heagerty P.J. and Zheng Y. (2005). Survival model predictive accuracy and ROC curves, Biometrics, 61:92–105.

Ishwaran H. and Kogalur U.B. (2007). Random survival forests for R, R News, 7(2):25–31.

Ishwaran H., Kogalur U.B., Blackstone E.H. and Lauer M.S. (2008). Random survival forests, Annals of Applied Statistics, 2:841–860.

See Also

plot.competing.risk.rfsrc, predict.rfsrc, rfsrc

Examples

## veteran data
data(veteran, package = "randomForestSRC")
plot.survival(rfsrc(Surv(time, status) ~ ., veteran),
              cens.model = "rfsrc")

## pbc data
data(pbc, package = "randomForestSRC")
pbc.obj <- rfsrc(Surv(days, status) ~ ., pbc)

## ------------------------------------------------------------
## Concordance error from the stored OOB risk predictions
## ------------------------------------------------------------
print(get.cindex(pbc.obj$yvar[, 1], pbc.obj$yvar[, 2],
                 pbc.obj$predicted.oob))

## standard survival diagnostics
plot.survival(pbc.obj)
plot.survival(pbc.obj, subset = c(3, 10), collapse = TRUE)

## Brier and AUCT helpers with pointwise intervals
brier.obj <- get.brier.survival(pbc.obj, conf.int = TRUE)
print(head(brier.obj$brier.score))
auct.obj <- get.auct.survival(pbc.obj, conf.int = TRUE)
print(head(auct.obj$auct.score))

## compare two grow-data censoring models
brier.km <- get.brier.survival(pbc.obj, cens.model = "km")
brier.rf <- get.brier.survival(pbc.obj, cens.model = "rfsrc")
plot(brier.km$brier.score$time,
     brier.km$brier.score$brier.score, type = "s", col = 2,
     xlab = "Time", ylab = "Brier Score")
lines(brier.rf$brier.score$time,
      brier.rf$brier.score$brier.score, type = "s", col = 4)
legend("bottomright",
       legend = c("cens.model = km", "cens.model = rfsrc"),
       col = c(2, 4), lty = 1)

## Brier and AUCT curves; Brier includes the grow-data KM null reference
perf <- plotBrierAUC(pbc.obj)
print(head(data.frame(
  time = perf$brier$time,
  forest = perf$brier$brier.score$brier.score,
  null = perf$brier$null$brier.score$brier.score
)))
plotBrierAUC(pbc.obj, plots = "auct", conf.int = 0.90)
plotBrierAUC(pbc.obj, plots = "brier", conf.int = 0.90,
             ylim = c(0, .4), null.lty = 3, legend.pos = "topleft")

## Omit the Brier reference, retaining the forest curve and its band
plotBrierAUC(pbc.obj, plots = "brier", brier.null = FALSE)

## Independent test evaluation: the null still uses training outcomes
set.seed(19)
pbc.complete <- na.omit(pbc)
trn <- sample(seq_len(nrow(pbc.complete)),
              size = floor(0.7 * nrow(pbc.complete)))
grow <- rfsrc(Surv(days, status) ~ ., pbc.complete[trn, ], ntree = 100)
test <- predict(grow, newdata = pbc.complete[-trn, ])
test.perf <- plotBrierAUC(test, plots = "brier")
print(c(forest = test.perf$brier$crps.std,
        null = test.perf$brier$null$crps.std))

## Obtain both curves without opening a graphics device
perf <- plotBrierAUC(test, plots = "brier", show.plots = FALSE)
print(head(perf$brier$null$brier.score))

randomForestSRC documentation built on Sept. 16, 2026, 5:06 p.m.