| Title: | Random Forest-Based Multistate Survival Analysis |
|---|---|
| Description: | Fits transition-specific cause-specific random survival forests on a clock-reset duration scale for acyclic, non-recurrent multistate processes. Entry-conditioned state-occupation probabilities are assembled from predicted cumulative hazards by semi-Markov entry-mass and sojourn convolution on a validated regular grid. The one-row-per-subject interface supports one common initial state, one recorded entry per state, baseline time-fixed covariates, competing exits, and independent right censoring. Left truncation, recurrent visits, directed cycles, time-dependent covariates, and ongoing-sojourn dynamic prediction are not supported. The package also provides calendar-time Aalen-Johansen point estimates as a covariate-free descriptive baseline, transition-specific permutation importance, genuine ranger edge OOB concordance, and patient-level cross-validated IPCW state-probability scoring. Methods are described in Ishwaran et al. (2008) <doi:10.1214/08-AOAS169> for random survival forests, Putter et al. (2007) <doi:10.1002/sim.2712> for multistate competing risks decomposition, and Aalen and Johansen (1978) <https://www.jstor.org/stable/4615704> for the nonparametric estimator. |
| Authors: | Yiqing Chen [aut, cre] |
| Maintainer: | Yiqing Chen <[email protected]> |
| License: | MIT + file LICENSE |
| Version: | 0.1.9 |
| Built: | 2026-10-11 15:05:04 UTC |
| Source: | https://github.com/evalunac/rfmstate |
Fits transition-specific cause-specific random survival forests on a clock-reset duration scale for acyclic, non-recurrent multistate processes. Patient/profile state probabilities are assembled by semi-Markov convolution and are conditional on fresh state entry. The package supports a common initial state, one recorded entry per state, baseline time-fixed covariates, competing exits, and independent right censoring. It does not support left truncation, cycles/recurrent visits, time-dependent covariates, or ongoing-sojourn dynamic prediction. A strict predictor contract prevents IDs, event/censoring times, response fields, and arbitrary long-format columns from entering a forest. Every full-data or cross-validation fit learns its predictor schema only from its own fitting rows, and every reported OOB statistic has verified OOB coverage. A calendar-time Aalen-Johansen point estimator is provided as a covariate-free descriptive baseline. The package provides:
State space and transition structure definition
Wide-to-long data conversion for multistate counting processes
Cause-specific random forest fitting per origin state
Entry-conditioned state probabilities via semi-Markov convolution
Aalen-Johansen point estimation (covariate-free baseline)
Per-transition feature importance
Genuine edge OOB concordance and patient-level cross-validated IPCW Brier scores
Comprehensive visualizations
Maintainer: Yiqing Chen [email protected]
Authors:
Yiqing Chen [email protected]
Useful links:
Computes nonparametric estimates of transition probabilities using the Aalen-Johansen estimator via Nelson-Aalen cumulative hazard increments and product-integral construction.
aalen_johansen(msdata, s = 0)aalen_johansen(msdata, s = 0)
msdata |
An |
s |
Legacy numeric starting-time argument. Only |
The Aalen-Johansen estimator generalizes the Kaplan-Meier estimator to multistate models. RFmstate exposes it as a descriptive, covariate-free calendar-time population benchmark from the recorded common baseline.
The estimator uses calendar-time risk sets and is separate from the package's clock-reset semi-Markov forest solver. It is a covariate-free descriptive benchmark from the common baseline. Observed competing exits are transitions, not external censoring. No clipping or row normalization is applied.
An object of class "aj_estimate" containing:
Numeric vector of unique event times.
List of calendar-time product-integral point-estimate matrices from the recorded common baseline at each event time.
Matrix of state occupation probabilities over time. Rows are time points, columns are states.
List of Nelson-Aalen cumulative hazard matrices.
List of hazard increment matrices at each event time.
Matrix of at-risk counts over time.
Data frame of event counts per transition.
The multistate structure used.
The starting time.
Calendar time from the common study origin.
Statement that no validated variance or confidence interval is returned.
Only the recorded common baseline s = 0 is supported. The function
does not support left truncation, general conditional for
s > 0, covariate adjustment, or validated variance/confidence bands.
The point estimate should not be described as the covariate-free form of the
clock-reset forest model.
ms <- clinical_states() dat <- sim_clinical_data(n = 200, structure = ms, seed = 42) msdata <- prepare_data( data = dat, id = "ID", structure = ms, time_map = list( Responded = "time_Responded", Unresponded = "time_Unresponded", Stabilized = "time_Stabilized", Progressed = "time_Progressed", Death = "time_Death" ), censor_col = "time_censored", covariates = c("age", "sex", "BMI", "treatment") ) aj <- aalen_johansen(msdata) print(aj)ms <- clinical_states() dat <- sim_clinical_data(n = 200, structure = ms, seed = 42) msdata <- prepare_data( data = dat, id = "ID", structure = ms, time_map = list( Responded = "time_Responded", Unresponded = "time_Unresponded", Stabilized = "time_Stabilized", Progressed = "time_Progressed", Death = "time_Death" ), censor_col = "time_censored", covariates = c("age", "sex", "BMI", "treatment") ) aj <- aalen_johansen(msdata) print(aj)
A convenience function that creates the standard clinical trial multistate structure with states: Baseline, Responded, Unresponded, Stabilized, Progressed, Death.
clinical_states()clinical_states()
An mstate_structure object.
This is a fixed demonstration DAG, not a claim that the package supports
cyclic or recurrent clinical histories. Use define_multistate() for
another supported single-root DAG.
ms <- clinical_states() print(ms)ms <- clinical_states() print(ms)
Combines one clock-reset cumulative cause-specific hazard curve per allowed
edge using a regular-grid entry-mass/sojourn convolution. The result is not
a general Markov matrix: each row is conditional on fresh entry
into its named starting state at elapsed duration zero.
compute_trans_prob( cum_hazards, structure, s = 0, times = NULL, start_state = NULL, grid_step = NULL, target_grid_points = 1024L, max_grid_points = 131073L, grid_tol = 5e-04, max_grid_refinements = 7L, check_grid = TRUE, prob_tol = 1e-08, extrapolate = c("error", "flat") )compute_trans_prob( cum_hazards, structure, s = 0, times = NULL, start_state = NULL, grid_step = NULL, target_grid_points = 1024L, max_grid_points = 131073L, grid_tol = 5e-04, max_grid_refinements = 7L, check_grid = TRUE, prob_tol = 1e-08, extrapolate = c("error", "flat") )
cum_hazards |
Named list of data frames, one for every allowed edge.
Names must be |
structure |
An |
s |
Legacy starting-time argument. Only zero is supported. |
times |
Finite nonnegative elapsed durations. If |
start_state |
State in which fresh entry occurs. The default is the
common initial state stored in |
grid_step |
Optional initial internal grid step. |
target_grid_points |
Initial number of regular-grid intervals when
|
max_grid_points |
Maximum number of internal grid points, including elapsed duration zero. |
grid_tol |
Maximum probability change allowed on grid refinement.
The default |
max_grid_refinements |
Maximum number of step halvings. |
check_grid |
Whether to require grid-refinement convergence. Keep
|
prob_tol |
Positive tolerance for probability bounds and unit mass. |
extrapolate |
Either |
Cumulative hazards are evaluated as right-continuous step functions. Within each regular interval, total exit probability is allocated to causes in proportion to their cumulative-hazard increments. Entry mass is assigned to the interval's right endpoint. Grid refinement controls these approximations. A positive cumulative hazard at duration zero is rejected because the data contract permits no instantaneous fresh-entry transition. Tiny monotonicity corrections within the declared tolerance are counted. Probabilities are never clipped or row-normalized.
A trans_prob object. entry_prob and its temporary
alias P are arrays ordered selected starting state by occupied
state by elapsed time; the starting-state dimension has length one.
state_occ is the selected starting-state slice. The
object also records time, structure, fresh-entry conditioning,
support and extrapolation policy, the final internal grid and step,
convergence error/refinements, roundoff corrections, and
prob_tol.
Curves must cover every declared edge and represent cumulative
cause-specific hazards on the clock-reset duration scale. The function does
not estimate hazards, accept left truncation or s > 0, condition on
an ongoing sojourn, support cyclic/recurrent graphs, or provide uncertainty.
Prediction beyond conservative support fails unless the explicit flat-hazard
sensitivity extension is selected.
ms <- define_multistate(c("A", "B"), "B", list(A = "B")) tt <- seq(0, 2, length.out = 2001) hazards <- list("A->B" = data.frame(time = tt, hazard = 0.4 * tt)) pr <- compute_trans_prob(hazards, ms, times = c(0, 1, 2), target_grid_points = 512) pr$state_occms <- define_multistate(c("A", "B"), "B", list(A = "B")) tt <- seq(0, 2, length.out = 2001) hazards <- list("A->B" = data.frame(time = tt, hazard = 0.4 * tt)) pr <- compute_trans_prob(hazards, ms, times = c(0, 1, 2), target_grid_points = 512) pr$state_occ
Defines the state space, absorbing states, and allowed transitions for a multistate model.
define_multistate(state_names, absorbing, transitions)define_multistate(state_names, absorbing, transitions)
state_names |
Character vector of state names. |
absorbing |
Character vector of absorbing state names (must be a subset
of |
transitions |
A named list where each element name is an origin state and the value is a character vector of destination states reachable from that origin. Absorbing states should not appear as list names. |
An object of class "mstate_structure" containing:
Character vector of all state names.
Integer, number of states.
Character vector of absorbing states.
Character vector of transient (non-absorbing) states.
Named list of allowed transitions.
Integer matrix where entry [i,j] is the
transition number for allowed transition i->j, or NA if not
allowed.
Total number of allowed transitions.
Data frame listing all transitions with columns
trans_id, from, to.
The unique graph root used as the common initial state.
State order used by the semi-Markov solver.
Only directed acyclic, non-recurrent structures with exactly one root and at least one absorbing state are supported. Every transient state must be reachable from the root and must have a directed path to absorption. Cycles, self-transitions, recurrent visits, disconnected states, and multiple subject-specific initial states are rejected.
ms <- define_multistate( state_names = c("Baseline", "Responded", "Progressed", "Death"), absorbing = "Death", transitions = list( Baseline = c("Responded", "Progressed", "Death"), Responded = c("Progressed", "Death"), Progressed = c("Death") ) ) print(ms)ms <- define_multistate( state_names = c("Baseline", "Responded", "Progressed", "Death"), absorbing = "Death", transitions = list( Baseline = c("Responded", "Progressed", "Death"), Responded = c("Progressed", "Death"), Progressed = c("Death") ) ) print(ms)
Reports genuine ranger edge-specific OOB concordance. Optionally, complete patient-level cross-validation refits every edge forest and returns IPCW full-state Brier scores and an integrated Brier score.
diagnose(object, ...) ## S3 method for class 'rfmstate' diagnose( object, eval_times = NULL, method = c("edge_oob", "cv"), folds = 5L, repeats = 1L, seed = 2026L, g_min = 0.05, ... )diagnose(object, ...) ## S3 method for class 'rfmstate' diagnose( object, eval_times = NULL, method = c("edge_oob", "cv"), folds = 5L, repeats = 1L, seed = 2026L, g_min = 0.05, ... )
object |
A fitted |
... |
Ignored; validation settings must use the documented arguments. |
eval_times |
Prespecified calendar-time validation horizons. For
|
method |
|
folds |
Positive integer number of patient-level folds for cross-validation. |
repeats |
Positive integer number of repeated fold assignments. |
seed |
Nonnegative integer fold-assignment and refit seed. |
g_min |
Minimum allowed training-fold censoring survival, strictly between zero and one. |
Ranger OOB concordance is 1 - prediction.error for each
binary cause-specific edge endpoint. It does not validate the assembled
multistate probability vector. Full-state scores use subject-level held-out
predictions and training-fold Kaplan–Meier censoring estimates. The
stronger marginal independent-censoring assumption applies to these KM
weights. Fold assignment occurs before schema construction; unordered
factor levels and numeric ranges are reconstructed from training subjects
only. A held-out-only factor level is an explicit fold failure. No
bias–variance decomposition is provided.
An rfmstate_diag object. Edge OOB tables
(edge_oob, oob_error, and concordance) are always
present. With cross-validation, brier is a full-state IPCW score
table, ibs is its trapezoidal integral over
integration_interval, and fold assignments/support/seeds are stored
in assignments and fold_summary; folds is a
backward-compatible alias for fold_summary. The object also records
evaluation times, whether their grid was prespecified or exploratory,
censoring model, g_min, method, and an explicit validation label.
Edge OOB concordance applies only to separate binary cause-specific
endpoints and is not a full-pipeline validation score. Cross-validation
currently uses a training-fold marginal Kaplan–Meier censoring model and
therefore requires marginal independent censoring for the score. Every fold
must fit every declared edge and cover all evaluation horizons; failures,
sparse edges, or censoring survival below g_min stop validation. No
calibration model, prediction interval, or bias–variance decomposition is
returned. Automatically generated evaluation times are exploratory and are
labeled as such; confirmatory work should supply prespecified times.
ms <- clinical_states() dat <- sim_clinical_data(300, structure = ms, seed = 42) long <- prepare_data( dat, "ID", ms, list(Responded="time_Responded", Unresponded="time_Unresponded", Stabilized="time_Stabilized", Progressed="time_Progressed", Death="time_Death"), "time_censored", c("age", "sex", "BMI", "treatment") ) fit <- rfmstate(long, num.trees = 100, min_events = 3, seed = 42) diagnose(fit)ms <- clinical_states() dat <- sim_clinical_data(300, structure = ms, seed = 42) long <- prepare_data( dat, "ID", ms, list(Responded="time_Responded", Unresponded="time_Unresponded", Stabilized="time_Stabilized", Progressed="time_Progressed", Death="time_Death"), "time_censored", c("age", "sex", "BMI", "treatment") ) fit <- rfmstate(long, num.trees = 100, min_events = 3, seed = 42) diagnose(fit)
Return ordinary data rows from an msdata object
## S3 method for class 'msdata' head(x, n = 6L, ...)## S3 method for class 'msdata' head(x, n = 6L, ...)
x |
An |
n |
Number of rows to return. |
... |
Ignored. |
An ordinary data.frame containing the first n rows.
Class-specific metadata are intentionally removed from the returned preview,
so it should not be passed to rfmstate(). The source msdata
object is unchanged.
ms <- define_multistate(c("A", "B"), "B", list(A = "B")) dat <- data.frame(id = 1:3, x = 1:3, time_B = 1:3, censor = NA_real_) long <- prepare_data(dat, "id", ms, list(B = "time_B"), "censor", "x") head(long, 2)ms <- define_multistate(c("A", "B"), "B", list(A = "B")) dat <- data.frame(id = 1:3, x = 1:3, time_B = 1:3, censor = NA_real_) long <- prepare_data(dat, "id", ms, list(B = "time_B"), "censor", "x") head(long, 2)
Extracts and organizes variable importance scores from the fitted random forest models for each transition. These are ranger permutation-importance scores based on edge-specific OOB predictive loss, not causal effects.
importance(object, ...) ## S3 method for class 'rfmstate' importance(object, ...)importance(object, ...) ## S3 method for class 'rfmstate' importance(object, ...)
object |
A fitted |
... |
Ignored. |
Values are transition-specific and may have different effective sample sizes. Negative values can arise from Monte Carlo variation, sparse events, correlated predictors, or noise; they do not indicate protective or causal effects. Repeated seeds/resamples are needed to assess stability.
An object of class "rfmstate_importance" containing:
Data frame with columns variable,
from, to, transition, n_events, and
importance.
Matrix with variables as rows and transitions as columns.
Covariate names.
Character vector of transition labels.
Importance describes transition-specific predictive contribution under the
fitted ranger endpoint. Values are not causal effects, are not directly
comparable across edges with different risk sets/event counts, and do not
validate the assembled full-state probability model. Models fitted with
importance = "none" return unavailable values.
ms <- clinical_states() dat <- sim_clinical_data(n = 200, structure = ms, seed = 42) msdata <- prepare_data( data = dat, id = "ID", structure = ms, time_map = list( Responded = "time_Responded", Unresponded = "time_Unresponded", Stabilized = "time_Stabilized", Progressed = "time_Progressed", Death = "time_Death" ), censor_col = "time_censored", covariates = c("age", "sex", "BMI", "treatment") ) fit <- rfmstate(msdata, num.trees = 100, seed = 42) imp <- importance(fit) print(imp)ms <- clinical_states() dat <- sim_clinical_data(n = 200, structure = ms, seed = 42) msdata <- prepare_data( data = dat, id = "ID", structure = ms, time_map = list( Responded = "time_Responded", Unresponded = "time_Unresponded", Stabilized = "time_Stabilized", Progressed = "time_Progressed", Death = "time_Death" ), censor_col = "time_censored", covariates = c("age", "sex", "BMI", "treatment") ) fit <- rfmstate(msdata, num.trees = 100, seed = 42) imp <- importance(fit) print(imp)
Draws a state transition diagram with event counts annotated on edges. Uses a layered layout that adapts to any number of states and automatically routes arrows around intermediate state boxes using Bezier curves when needed.
plot_transition_diagram( structure, msdata = NULL, col = NULL, main = "Transition Diagram", xlab = "", ylab = "", ... )plot_transition_diagram( structure, msdata = NULL, col = NULL, main = "Transition Diagram", xlab = "", ylab = "", ... )
structure |
An |
msdata |
Optional |
col |
Node colors. Default uses the standard palette. |
main |
Title. |
xlab, ylab
|
Optional axis labels; defaults are blank because graph coordinates have no statistical scale. |
... |
Ignored. |
No return value, called for its side effect of producing a plot.
The layout is descriptive and preserves the validated DAG. It does not
represent cycles/recurrent visits, transition intensity, uncertainty, or
edge direction over calendar time. Counts are shown only when a compatible
msdata object is supplied.
ms <- clinical_states() plot_transition_diagram(ms, xlab = "Graph layer")ms <- clinical_states() plot_transition_diagram(ms, xlab = "Graph layer")
Visualizes calendar-time state occupation, product-integral point estimates, Nelson–Aalen cumulative cause-specific hazards, or hazard increments from the Aalen–Johansen estimator.
## S3 method for class 'aj_estimate' plot( x, type = c("state_occupation", "stacked_transition_prob", "cumulative_hazard", "hazard_increment", "transition_intensity"), states = NULL, ci = FALSE, col = NULL, main = NULL, xlab = "Time", ylab = NULL, ... )## S3 method for class 'aj_estimate' plot( x, type = c("state_occupation", "stacked_transition_prob", "cumulative_hazard", "hazard_increment", "transition_intensity"), states = NULL, ci = FALSE, col = NULL, main = NULL, xlab = "Time", ylab = NULL, ... )
x |
An |
type |
Character, one of |
states |
Character vector of states to plot (default: all). For cumulative-hazard and hazard-increment plots, transitions are filtered by destination state. |
ci |
One nonmissing logical value. This is a deprecated compatibility
argument: |
col |
One or more colors for selected states/transitions. Values are
recycled when necessary; |
main |
Title (default: auto-generated). |
xlab, ylab
|
Axis labels. |
... |
Additional arguments passed to |
The input x object, returned invisibly. Called for its
side effect of producing a plot.
AJ plots are covariate-free calendar-time point estimates from the common
baseline. Confidence bands are unavailable. Hazard increments are discrete
Nelson–Aalen increments, not smoothed transition intensities. The legacy
"transition_intensity" type is deprecated and draws the same
hazard-increment plot with a warning. The historical
"stacked_transition_prob" label displays the selected state-occupation
components from the recorded initial state.
ms <- define_multistate(c("A", "B"), "B", list(A = "B")) dat <- data.frame(id = 1:6, x = 1:6, time_B = 1:6, censor = NA_real_) long <- prepare_data(dat, "id", ms, list(B = "time_B"), "censor", "x") aj <- aalen_johansen(long) plot(aj, states = "B", xlab = "Study day", ylab = "Probability") plot(aj, type = "hazard_increment", states = "B")ms <- define_multistate(c("A", "B"), "B", list(A = "B")) dat <- data.frame(id = 1:6, x = 1:6, time_B = 1:6, censor = NA_real_) long <- prepare_data(dat, "id", ms, list(B = "time_B"), "censor", "x") aj <- aalen_johansen(long) plot(aj, states = "B", xlab = "Study day", ylab = "Probability") plot(aj, type = "hazard_increment", states = "B")
Visualizes genuine edge OOB concordance or patient-level cross-validated full-state Brier scores.
## S3 method for class 'rfmstate_diag' plot( x, type = c("brier", "concordance"), col = NULL, main = NULL, xlab = NULL, ylab = NULL, ... )## S3 method for class 'rfmstate_diag' plot( x, type = c("brier", "concordance"), col = NULL, main = NULL, xlab = NULL, ylab = NULL, ... )
x |
An |
type |
Character, one of |
col |
One or more plot colors; recycled when necessary. |
main |
Title. |
xlab, ylab
|
Axis labels. |
... |
Additional arguments. |
The input x object, returned invisibly. Called for its
side effect of producing a plot.
Concordance plots contain separate edge-level ranger OOB statistics, not
assembled full-state validation. Brier plots require an object produced by
diagnose(..., method = "cv") and inherit its censoring/support
limitations. No bias–variance plot is available.
d <- structure(list( concordance = data.frame(transition = "A->B", c_index = 0.7), brier = NULL ), class = "rfmstate_diag") plot(d, type = "concordance", xlab = "Edge", ylab = "OOB C-index")d <- structure(list( concordance = data.frame(transition = "A->B", c_index = 0.7), brier = NULL ), class = "rfmstate_diag") plot(d, type = "concordance", xlab = "Edge", ylab = "OOB C-index")
Visualizes per-transition feature importance as a grouped barplot or heatmap.
## S3 method for class 'rfmstate_importance' plot( x, type = c("barplot", "heatmap"), col = NULL, main = NULL, xlab = NULL, ylab = NULL, ... )## S3 method for class 'rfmstate_importance' plot( x, type = c("barplot", "heatmap"), col = NULL, main = NULL, xlab = NULL, ylab = NULL, ... )
x |
An |
type |
Character, one of |
col |
Bar colors or heatmap gradient-anchor colors. |
main |
Title. |
xlab, ylab
|
Axis labels. |
... |
Additional graphical arguments. For a heatmap they are passed to
|
The input x object, returned invisibly. Called for its
side effect of producing a plot.
Values are transition-specific ranger importance scores. Negative values may arise from Monte Carlo noise, sparse events, correlated predictors, or irrelevant variables and are not protective or causal effects. Scales may differ across transitions.
imp <- structure(list( importance_matrix = matrix( c(0.2, -0.1), nrow = 2, dimnames = list(c("x", "z"), "A->B") ) ), class = "rfmstate_importance") plot(imp, xlab = "OOB loss increase", ylab = "Predictor")imp <- structure(list( importance_matrix = matrix( c(0.2, -0.1), nrow = 2, dimnames = list(c("x", "z"), "A->B") ) ), class = "rfmstate_importance") plot(imp, xlab = "OOB loss increase", ylab = "Predictor")
Visualizes predicted state occupation probabilities and transition probabilities for individual patients.
## S3 method for class 'rfmstate_pred' plot( x, type = c("state_occupation", "transition_prob"), subject = 1L, col = NULL, main = NULL, states = NULL, xlab = "Elapsed duration", ylab = NULL, ... )## S3 method for class 'rfmstate_pred' plot( x, type = c("state_occupation", "transition_prob"), subject = 1L, col = NULL, main = NULL, states = NULL, xlab = "Elapsed duration", ylab = NULL, ... )
x |
An |
type |
Character, one of |
subject |
Integer, which subject to plot (default 1). Use 0 for mean across all subjects. |
col |
One or more state colors. Values are recycled when necessary;
|
main |
Title. |
states |
Occupied states to display. |
xlab, ylab
|
Axis labels. |
... |
Additional arguments passed to |
The input x object, returned invisibly. Called for its
side effect of producing a plot.
Curves are entry-conditioned on fresh entry into the fitted
start_state; they are not general Markov curves and have
no confidence bands. subject = 0 is an arithmetic mean of profile
predictions, not a population-standardized estimator. The historical
"transition_prob" plot type displays the occupied-state components
conditional on the selected fresh-entry state.
ms <- define_multistate(c("A", "B"), "B", list(A = "B")) dat <- data.frame(id = 1:60, x = seq(-1, 1, length.out = 60), time_B = 1:60, censor = NA_real_) long <- prepare_data(dat, "id", ms, list(B = "time_B"), "censor", "x") fit <- rfmstate(long, num.trees = 20, min.node.size = 3, min_events = 3, sparse_warning = Inf, seed = 42) pred <- predict(fit, data.frame(x = 0), times = c(0, 5, 10), target_grid_points = 64, max_grid_points = 1025) plot(pred, states = c("A", "B"), xlab = "Elapsed day") plot(pred, type = "transition_prob", states = "B")ms <- define_multistate(c("A", "B"), "B", list(A = "B")) dat <- data.frame(id = 1:60, x = seq(-1, 1, length.out = 60), time_B = 1:60, censor = NA_real_) long <- prepare_data(dat, "id", ms, list(B = "time_B"), "censor", "x") fit <- rfmstate(long, num.trees = 20, min.node.size = 3, min_events = 3, sparse_warning = Inf, seed = 42) pred <- predict(fit, data.frame(x = 0), times = c(0, 5, 10), target_grid_points = 64, max_grid_points = 1025) plot(pred, states = c("A", "B"), xlab = "Elapsed day") plot(pred, type = "transition_prob", states = "B")
Predicts clock-reset semi-Markov state-occupation probabilities conditional on fresh entry into a selected state at elapsed duration zero.
## S3 method for class 'rfmstate' predict( object, newdata = NULL, times = NULL, s = 0, start_state = NULL, grid_step = NULL, target_grid_points = 1024L, max_grid_points = 131073L, grid_tol = 5e-04, max_grid_refinements = 7L, check_grid = TRUE, extrapolate = c("error", "flat"), ... )## S3 method for class 'rfmstate' predict( object, newdata = NULL, times = NULL, s = 0, start_state = NULL, grid_step = NULL, target_grid_points = 1024L, max_grid_points = 131073L, grid_tol = 5e-04, max_grid_refinements = 7L, check_grid = TRUE, extrapolate = c("error", "flat"), ... )
object |
A fitted |
newdata |
Prediction profiles. |
times |
Finite nonnegative elapsed durations. |
s |
Legacy start-time argument; only zero is supported. |
start_state |
Fresh-entry starting state. Defaults to the common initial state. |
grid_step |
Optional positive initial internal grid step. |
target_grid_points |
Initial regular-grid interval count when
|
max_grid_points |
Maximum allowed internal grid-point count. |
grid_tol |
Maximum change permitted between successive refinements. |
max_grid_refinements |
Maximum number of grid-step halvings. |
check_grid |
Whether convergence is required. Keep |
extrapolate |
|
... |
Ignored; no ranger arguments are accepted during prediction. |
These are not general Markov matrices. Predictions do
not condition on an already elapsed sojourn and carry no confidence bands.
By default the horizon cannot exceed the minimum observed sojourn support
among reachable transient states.
An rfmstate_pred object. entry_prob and the temporary
alias P have dimensions profile by selected starting state by
occupied state by elapsed time; the starting-state dimension has length
one. state_occ is that selected starting-state slice.
The object also contains requested time, per-profile edge
cum_hazard, validated newdata, used columns, support and
extrapolation metadata, grid-convergence records, conditioning statement,
predictor-support extrapolation records, fit-specific schema identifier,
state structure, and package versions.
Prediction is conditional on fresh entry into start_state at elapsed
duration zero. General prediction for s > 0, left
truncation, ongoing-sojourn/landmark prediction, recurrent histories, and
confidence intervals are unavailable. Missing covariates, nonfinite numeric
values, incompatible classes, and unseen factor levels are rejected.
Numeric values outside the fit-specific range warn and are recorded as
predictor-support extrapolation; they are never truncated.
Times beyond support fail unless extrapolate = "flat"; that option
assumes zero additional hazard and is a sensitivity analysis only.
ms <- clinical_states() dat <- sim_clinical_data(300, structure = ms, seed = 42) long <- prepare_data( dat, "ID", ms, list(Responded="time_Responded", Unresponded="time_Unresponded", Stabilized="time_Stabilized", Progressed="time_Progressed", Death="time_Death"), "time_censored", c("age", "sex", "BMI", "treatment") ) fit <- rfmstate(long, num.trees = 100, min_events = 3, seed = 42) horizon <- min(fit$max_duration_by_origin) predict(fit, data.frame(age=60, sex=1, BMI=26, treatment=1), times = c(0, horizon / 2))ms <- clinical_states() dat <- sim_clinical_data(300, structure = ms, seed = 42) long <- prepare_data( dat, "ID", ms, list(Responded="time_Responded", Unresponded="time_Unresponded", Stabilized="time_Stabilized", Progressed="time_Progressed", Death="time_Death"), "time_censored", c("age", "sex", "BMI", "treatment") ) fit <- rfmstate(long, num.trees = 100, min_events = 3, seed = 42) horizon <- min(fit$max_duration_by_origin) predict(fit, data.frame(age=60, sex=1, BMI=26, treatment=1), times = c(0, horizon / 2))
Converts one-row-per-subject wide data into validated, long-form sojourn records for an acyclic, non-recurrent multistate process. Entry times are calendar times from a common origin; forest response times are the resulting clock-reset durations within states.
prepare_data( data, id, structure, time_map, censor_col, covariates, initial_state = NULL )prepare_data( data, id, structure, time_map, censor_col, covariates, initial_state = NULL )
data |
A data frame with one row per subject. |
id |
Name of the unique subject-ID column. |
structure |
An |
time_map |
Named list mapping every noninitial state exactly once to its first and only calendar-time entry column. The initial state is excluded. |
censor_col |
Name of the external right-censoring column. |
covariates |
Explicit character vector of baseline, time-fixed
predictors. Use |
initial_state |
Common initial state. It must equal the unique graph
root; the default is stored in |
Delayed entry, recurrent visits, tied entry times, time-dependent
covariates, and missing covariates are unsupported. Invalid trajectories
produce an error; no event or interval is silently discarded. Character
covariates are converted to factors in first-observed level order, and the
resulting preparation schema is descriptive only. Every model fit and
cross-validation refit rebuilds its schema from its own fitting rows.
Subject ID, mapped entry times, censoring time, reserved long-format names,
and names beginning .rfm_ cannot be predictors.
An msdata data frame containing id, from,
to, Tstart, Tstop, status,
trans_id, duration, and retained covariates. Metadata store
the graph, source-role table, predictor contract, descriptive preparation
schema, counts, and observed duration support.
The input must contain one row per subject, one common initial state at time zero, and at most one exact entry time per noninitial state. Left truncation, cycles/recurrent visits, tied or interval-censored transitions, time-dependent covariates, missing fitting covariates, events after censoring, and events after absorption are rejected. A nonabsorbed subject requires a finite censoring time strictly after the last state entry.
ms <- clinical_states() dat <- sim_clinical_data(n = 50, structure = ms, seed = 42) msdata <- prepare_data( dat, id = "ID", structure = ms, time_map = list( Responded = "time_Responded", Unresponded = "time_Unresponded", Stabilized = "time_Stabilized", Progressed = "time_Progressed", Death = "time_Death" ), censor_col = "time_censored", covariates = c("age", "sex", "BMI", "treatment") ) head(msdata)ms <- clinical_states() dat <- sim_clinical_data(n = 50, structure = ms, seed = 42) msdata <- prepare_data( dat, id = "ID", structure = ms, time_map = list( Responded = "time_Responded", Unresponded = "time_Unresponded", Stabilized = "time_Stabilized", Progressed = "time_Progressed", Death = "time_Death" ), censor_col = "time_censored", covariates = c("age", "sex", "BMI", "treatment") ) head(msdata)
Concise print methods for fitted models, predictions, direct semi-Markov probability results, Aalen–Johansen point estimates, importance results, diagnostics, and fitted-model summaries.
## S3 method for class 'aj_estimate' print(x, ...) ## S3 method for class 'rfmstate_diag' print(x, ...) ## S3 method for class 'rfmstate' print(x, ...) ## S3 method for class 'rfmstate_importance' print(x, ...) ## S3 method for class 'rfmstate_pred' print(x, ...) ## S3 method for class 'summary.rfmstate' print(x, ...) ## S3 method for class 'trans_prob' print(x, ...)## S3 method for class 'aj_estimate' print(x, ...) ## S3 method for class 'rfmstate_diag' print(x, ...) ## S3 method for class 'rfmstate' print(x, ...) ## S3 method for class 'rfmstate_importance' print(x, ...) ## S3 method for class 'rfmstate_pred' print(x, ...) ## S3 method for class 'summary.rfmstate' print(x, ...) ## S3 method for class 'trans_prob' print(x, ...)
x |
An object of the class required by the selected print method:
|
... |
Ignored. |
Printed values are deliberately concise and rounded only for
display. Stored numerical results are unchanged. The fitted-model and
summary methods identify edge-level ranger OOB quantities; prediction and
direct-probability methods state their fresh-entry conditioning; the AJ
method states that it returns point estimates only; and the diagnostic
method distinguishes edge OOB output from patient-level cross-validation.
The example is in a donttest block because it fits survival forests to
exercise every result class.
x, invisibly.
Print methods are descriptive and do not revalidate, refit, or add uncertainty. Rounded importance or probability values must not be used in downstream calculations. Edge OOB statistics are not full-state validation, negative permutation importance is not a protective/causal effect, forest confidence intervals are unavailable, and AJ confidence intervals are not returned.
ms <- define_multistate(c("A", "B"), "B", list(A = "B")) dat <- data.frame(id = 1:60, x = seq(-1, 1, length.out = 60), time_B = 1:60, censor = NA_real_) long <- prepare_data(dat, "id", ms, list(B = "time_B"), "censor", "x") fit <- rfmstate(long, num.trees = 20, min.node.size = 3, min_events = 3, sparse_warning = Inf, seed = 42) pred <- predict(fit, data.frame(x = 0), times = c(0, 10), target_grid_points = 64, max_grid_points = 1025) tt <- seq(0, 10, length.out = 1001) tp <- compute_trans_prob( list("A->B" = data.frame(time = tt, hazard = 0.1 * tt)), ms, times = c(0, 5, 10), target_grid_points = 128 ) aj <- aalen_johansen(long) print(fit) print(pred) print(tp) print(aj) print(importance(fit)) print(diagnose(fit)) print(summary(fit))ms <- define_multistate(c("A", "B"), "B", list(A = "B")) dat <- data.frame(id = 1:60, x = seq(-1, 1, length.out = 60), time_B = 1:60, censor = NA_real_) long <- prepare_data(dat, "id", ms, list(B = "time_B"), "censor", "x") fit <- rfmstate(long, num.trees = 20, min.node.size = 3, min_events = 3, sparse_warning = Inf, seed = 42) pred <- predict(fit, data.frame(x = 0), times = c(0, 10), target_grid_points = 64, max_grid_points = 1025) tt <- seq(0, 10, length.out = 1001) tp <- compute_trans_prob( list("A->B" = data.frame(time = tt, hazard = 0.1 * tt)), ms, times = c(0, 5, 10), target_grid_points = 128 ) aj <- aalen_johansen(long) print(fit) print(pred) print(tp) print(aj) print(importance(fit)) print(diagnose(fit)) print(summary(fit))
Prints a concise validation summary for an msdata object, including
patient and interval totals and an observed-transition table in the original
state-definition order. Use head() to display actual data rows.
## S3 method for class 'msdata' print(x, ...)## S3 method for class 'msdata' print(x, ...)
x |
An |
... |
Ignored. |
x, invisibly.
The printout is descriptive and does not refit or revalidate a modified
object. Re-run prepare_data() after changing rows or metadata.
ms <- define_multistate(c("A", "B"), "B", list(A = "B")) dat <- data.frame(id = 1:3, x = 1:3, time_B = 1:3, censor = NA_real_) long <- prepare_data(dat, "id", ms, list(B = "time_B"), "censor", "x") print(long)ms <- define_multistate(c("A", "B"), "B", list(A = "B")) dat <- data.frame(id = 1:3, x = 1:3, time_B = 1:3, censor = NA_real_) long <- prepare_data(dat, "id", ms, list(B = "time_B"), "censor", "x") print(long)
Prints the user-defined display order, absorbing states, common initial state, computational topological order, and numbered allowed transitions.
## S3 method for class 'mstate_structure' print(x, ...)## S3 method for class 'mstate_structure' print(x, ...)
x |
An |
... |
Ignored. |
x, invisibly.
The printout summarizes the validated graph; it does not imply support for cycles, recurrent visits, left truncation, or time-dependent covariates.
ms <- define_multistate(c("A", "B"), "B", list(A = "B")) print(ms)ms <- define_multistate(c("A", "B"), "B", list(A = "B")) print(ms)
Fits one ranger cause-specific survival forest for every declared edge of a validated acyclic, non-recurrent multistate structure. Forest response time is duration since fresh entry into the origin state.
rfmstate( msdata, covariates = NULL, num.trees = 1000L, mtry = NULL, min.node.size = 15L, min_events = 5L, sparse_warning = 50L, importance = "permutation", seed = NULL, ... )rfmstate( msdata, covariates = NULL, num.trees = 1000L, mtry = NULL, min.node.size = 15L, min_events = 5L, sparse_warning = 50L, importance = "permutation", seed = NULL, ... )
msdata |
An |
covariates |
Predictor names. |
num.trees |
Positive integer number of trees fitted for every edge. |
mtry |
Positive integer number of predictors considered at each split.
|
min.node.size |
Positive integer ranger minimum terminal-node size. |
min_events |
Explicit technical safeguard for target events per edge. It is not a universal adequacy threshold. |
sparse_warning |
Positive integer descriptive event-count threshold for
an edge warning; set |
importance |
One of |
seed |
|
... |
Supported ranger controls: |
Every competing exit from an origin state remains an observed exit
for that sojourn, but is coded as a non-target outcome in the binary
cause-specific forest for a particular edge. An unestimable declared edge
stops the whole fit; it is never omitted or represented by zero hazard.
The ranger model frame is constructed anew as .rfm_time,
.rfm_event, and the contract-approved predictors. Predictor factor
levels and numeric ranges are learned from the actual fitting rows, not
copied from a preparation-time or full-data schema.
An rfmstate object containing models for every
declared edge; structure; selected covariates and
predictor_schema; origin-state fitting data; backend event-time
grids; per-edge event, support, ranger-argument, and genuine OOB metadata;
the validated msdata; fit params; package versions; verified
OOB coverage; and the clock-reset semi-Markov time-scale contract.
RFmstate controls and rejects duplicate specification of formula,
data, num.trees, mtry, min.node.size,
importance, seed, write.forest, oob.error, and
keep.inbag.
Unnamed and unknown arguments also fail. The only names accepted through
... are:
Bootstrap/subsampling controls.
Survival split controls supported by the installed ranger version.
Unordered-factor handling.
Computation controls.
Maximum tree depth.
Predictors always considered for splitting.
These arguments retain ranger's definitions and are checked against the
installed ranger formal arguments before fitting. case.weights,
class.weights, split.select.weights, response controls, and
every other unlisted ranger argument are rejected in this release.
Effective sampling defaults are resolved and stored. In particular,
replace = FALSE, sample.fraction = 1 is rejected because it leaves no
OOB observations. RFmstate forces oob.error = TRUE and
keep.inbag = TRUE, verifies finite ranger OOB error after every edge
fit, and stores per-sojourn OOB-tree coverage.
The fit supports baseline, time-fixed, complete covariates in single-root
acyclic non-recurrent data. It does not support left truncation,
time-dependent covariates, recurrent visits/cycles, missing fitting
covariates, subject-specific frailty, clock-forward hazards, or confidence
intervals. Structural, outcome-time, censoring, ID, and arbitrary
long-format columns cannot be added as predictors. Every declared edge must
meet min_events; this technical
safeguard is not a universal adequacy threshold.
ms <- clinical_states() dat <- sim_clinical_data(300, structure = ms, seed = 42) long <- prepare_data( dat, "ID", ms, list(Responded="time_Responded", Unresponded="time_Unresponded", Stabilized="time_Stabilized", Progressed="time_Progressed", Death="time_Death"), "time_censored", c("age", "sex", "BMI", "treatment") ) fit <- rfmstate(long, num.trees = 100, min_events = 3, seed = 42)ms <- clinical_states() dat <- sim_clinical_data(300, structure = ms, seed = 42) long <- prepare_data( dat, "ID", ms, list(Responded="time_Responded", Unresponded="time_Unresponded", Stabilized="time_Stabilized", Progressed="time_Progressed", Death="time_Death"), "time_censored", c("age", "sex", "BMI", "treatment") ) fit <- rfmstate(long, num.trees = 100, min_events = 3, seed = 42)
Generates realistic clinical trial data with covariates and multistate event times for testing and demonstration. The structure must be an acyclic, non-recurrent graph with one common initial state.
sim_clinical_data(n = 500, structure = NULL, max_followup = 365, seed = NULL)sim_clinical_data(n = 500, structure = NULL, max_followup = 365, seed = NULL)
n |
Integer, number of patients to simulate. |
structure |
An |
max_followup |
Numeric, maximum follow-up time (for generating censoring). Default 365. |
seed |
Optional integer for reproducibility. |
Cause-specific transition hazards follow Weibull distributions with covariate effects
on the scale parameter. For the default clinical_states()
structure, transition-specific parameters are calibrated to produce
realistic clinical trial trajectories. For custom structures, sensible
default parameters are used for all transitions.
Event waiting times remain at full numerical precision. External
censoring is generated before each path and truncates it, so no event is
retained after censoring or after entry into any absorbing state. Supplying
the same seed, arguments, structure, and package version reproduces
the returned data.
A data frame in wide format with columns:
Patient identifier (1 to n).
Continuous, simulated from Normal(60, 12).
Binary 0/1.
Continuous, simulated from Normal(26, 5).
Binary 0/1 (balanced arms).
For each non-initial state in the
structure, the time (days) of entry into that state, or NA
if the state was not visited. Column names follow the pattern
time_<StateName> (e.g., time_Death).
Days until last follow-up (right censoring
time), or NA if an absorbing state was reached.
This helper supplies demonstration and test data, not the final manuscript simulation design. It uses four fixed baseline covariates and simple transition-specific Weibull cause-specific hazards. It does not generate left truncation, recurrent/cyclic histories, time-dependent covariates, interval-censored transitions, or ongoing-sojourn predictions.
dat <- sim_clinical_data(n = 100, seed = 123) head(dat) summary(dat)dat <- sim_clinical_data(n = 100, seed = 123) head(dat) summary(dat)
Provides a comprehensive summary of the fitted model including per-origin state information, OOB prediction error, and transition event counts.
## S3 method for class 'rfmstate' summary(object, ...)## S3 method for class 'rfmstate' summary(object, ...)
object |
A fitted |
... |
Ignored. |
The per-edge OOB concordance is 1 - prediction.error from
ranger for that edge's binary cause-specific endpoint. The printed summary
also lists every forwarded ranger argument.
An object of class "summary.rfmstate" containing the stored
call, covariates, forest parameters, patient/event/interval totals,
structure, and a per-edge table of risk-set size, target events,
competing/non-target exits, external censoring, ranger OOB error,
ranger OOB concordance, and verified OOB coverage.
Edge OOB values do not validate assembled full-state probabilities and may
be unstable for sparse transitions. Use diagnose(..., method = "cv")
for patient-level held-out full-state Brier scoring.
ms <- clinical_states() dat <- sim_clinical_data(n = 200, structure = ms, seed = 42) msdata <- prepare_data( data = dat, id = "ID", structure = ms, time_map = list( Responded = "time_Responded", Unresponded = "time_Unresponded", Stabilized = "time_Stabilized", Progressed = "time_Progressed", Death = "time_Death" ), censor_col = "time_censored", covariates = c("age", "sex", "BMI", "treatment") ) fit <- rfmstate(msdata, num.trees = 100, seed = 42) summary(fit)ms <- clinical_states() dat <- sim_clinical_data(n = 200, structure = ms, seed = 42) msdata <- prepare_data( data = dat, id = "ID", structure = ms, time_map = list( Responded = "time_Responded", Unresponded = "time_Unresponded", Stabilized = "time_Stabilized", Progressed = "time_Progressed", Death = "time_Death" ), censor_col = "time_censored", covariates = c("age", "sex", "BMI", "treatment") ) fit <- rfmstate(msdata, num.trees = 100, seed = 42) summary(fit)