3 Event Prediction: event_prediction
3.1 Introduction
The event_prediction function uses a combination of fit_tte_data and the numerical integration approaches found in nph_traj to produce estimates for event numbers at any given time point. There are four main inputs required for event prediction: Firstly a patient-level data set from which an event Curve can be estimated (‘data’). Secondly, an RCurve reflecting the combined already-observed recruitment and predicted future recruitment (‘rcurve’). Thirdly, a Curve reflecting the dropout censoring in the trial (‘dcurve’). Finally, optionally, a set of parameters describing what has actually been observed at an observed time point, to allow for conditional predictions; these are described in section 3.5.
Particular care should be paid to time units since measured event times are typically in days, while recruitment and timeline planning are typically in months (see section 3.2 for more on time units). event_prediction assumes that all times are in units of months, except for the event times, which may be in units of days (default) or months; this choice is controlled by the optional ‘units’ argument.
3.2 Event Curve Fitting
Section 2.3 describes the MLE fitting process and the arguments required for this. All of these arguments carry over into event_prediction.
3.3 Recruitment Specification
The RCurve for an event prediction should have units of months and typically be of type PieceR since the exact numbers of patients that have so far entered the trial in each month will be known. The program does not consider variability in future recruitment as it is expected that by the point that event prediction is performed a sizeable amount of recruitment will have already occurred, and the variability in future recruitment will be both low and subject to negative feedback (i.e. if it drops, efforts are made to increase it and vice versa). As event prediction is performed in a blinded manner, the randomisation ratio (or even the number of arms) is irrelevant, and therefore any between-arm distribution of patients can be specified.
3.4 Censoring Specification
The amount of dropout censoring varies wildly between trials. The endpoint Overall Survival (OS) ought to have minimal dropout that does not tangibly affect event prediction, whereas other endpoints are likely to have enough that it warrants consideration. A dropout distribution can be supplied to event_prediction via the ‘dcurve’ argument.
It is suggested to use fit_KM or fit_tte_data to curve-fit a version of the data where dropouts are handled as events, while administrative censoring and events are handled as censored. This should produce a reasonable estimate for censoring unconditional on event occurrence (assuming independent censoring).
3.5 Conditioning
To calculate conditional predictions, it is necessary to specify the number of observed events at a given time, as well as preferably the number of patients remaining at risk (total number intended to be recruited minus events occurred, minus dropouts). The conditional number at risk is important because if it is known that a larger-than-expected number of patients have dropped out, then it is clear this will lower the future expected rate of events (as fewer patients than expected are at risk). Conditional predictions can be enabled by specifying a conditional time through the ‘cond_Time’argument. The number of conditional events is then given by ’cond_Events’ and the conditional number at risk by ‘cond_NatRisk’. If the latter is not known, gestate calculates an estimate based upon the expected number at risk at the conditioning time adjusted by the difference between the conditional and expected numbers of events at that time.
3.6 Accounting for Effects of Adjudication
For trials with adjudication it is common for event prediction to be performed upon investigator-reported data due to delays in adjudication. However, typically a proportion of adjudicated events will be ‘down-adjudicated’, resulting in events used in the prediction not being ‘valid’. To account for this net down-adjudication (and the rarer case of net up-adjudication), the ‘discountHR’ argument is available. If it is assumed that patients remain at risk after an observed event, ‘discountHR’ corresponds to the probability of an event being adjudicated as true. It applies a discount hazard ratio to the fitted event curve to reflect this. The effect of specifying a value below 1 for this parameter is to ensure there are fewer events predicted initially, but since the patients remain at risk, the predictions after the discount is applied will slowly catch up with those without the discount. This approach can also be used in cases of up-adjudication where patients may be found to retrospectively have had an event at a later stage - here the value should be above 1.
The argument’discountHR’ may only be used in conjunction with a Weibull curve since Lognormal curves are not subject to proportional hazards. If patients are no longer at risk following an event (e.g. if trial participation is ended) but they may still be down-adjudicated, then this may instead be handled by simply multiplying all predicted event numbers by the probability of adjudication as an event.
If the ‘discountHR’ argument is used, care should be taken with conditional predictions; in this situation, the number of events to be conditioned on is the unadjusted (i.e. ‘observed’) number of events, rather than the adjusted number. Likewise, the conditioned number at risk should be for the unadjusted events rather than the adjusted. The program will automatically adjust the conditioning values to reflect the adjusted scale before performing the calculations.
3.7 Calculations
A Curve object is created from the fitted parameters. This is combined with the recruitment RCurve and censoring Curve to produce estimates of the expected event numbers at each integer time point using the same approach as nph_traj. Conditioned trajectories are created by scaling future events by the ratio of the conditioned number at risk to the expected number at risk at the conditioning time. The scaled number of additional events is then added to the conditioned number of events. Where no conditioned number at risk is provided, gestate assumes by default that it is the expected number minus the difference between conditioned and expected event numbers at the conditioning time; i.e. that each additional observed event than expected reduces the expected number at risk by one.
The estimation errors for the parameters are also carried over into the event prediction calculation and propagated through by the delta method to enable the estimation of standard errors on both the conditional and unconditional expected trajectories. The conditional standard error is always lower at any given time, particularly close to the time of conditioning (and it is 0 at that time). The standard errors of the trajectories are then combined with the number at risk and the event numbers via a beta binomial distribution to create prediction intervals for the conditional and unconditional predicted trajectories.
The ‘PI’ argument specifies the width of prediction interval to produce, by default 0.95 (95%) intervals.
A paper on the general approach to analytic prediction intervals and the formulae for the calculations is to be published.
3.8 Plotting Event Prediction
For convenience, there is an inbuilt plotting function for event prediction output; plot_ep. This is designed to automatically read output from event_prediction, so only requires the ‘data’ argument, specifying the name of the event_prediction output. By default, it makes use of the available information, plotting the unconditional trajectory with prediction CI, as well as the conditional ones if available. The arguments ‘trajectory’ and ‘which_PI’ can be used to eliminate some or all of these lines in the plot. The prediction_fitting argument switches the prediction intervals for confidence intervals of fitting for the mean trajectory estimate, which can be used for technical/diagnostic purposes; these parameter fitting CIs are not representative of the range of events that could occur.
To plot the observed events as well, an observed event trajectory may also be supplied through the ‘observed’ argument. This accepts either a vector of events, corresponding to the values at times 1,2,3…, or a 2-column matrix / data-frame with the first column the time points in ascending order, and the second the event numbers.
A target number of events may also be specified using the ‘target’ argument. The optional ‘max_time’ and ‘max_E’ options override the automatic plotting area dimensions, while the ‘legend_position’ and ‘no_legend’ arguments control the legend positioning and whether it is displayed at all. Standard graphical arguments are also accepted.
3.9 Diagnostics: Plotting the Fitted Kaplan Meier Curve
To assess the goodness of fit for the event fitting, it is often helpful to compare the fitted event distribution to the observed Kaplan Meier plot. The function plot_km_fit will automatically create such a plot from the relevant input. It takes two key arguments: ‘fit’ and ‘data’. ‘fit’ should be the name of the output from fit_KM, fit_tte_data, or event_prediction.’data’ should be the name of the original patient data. The arguments ‘Time’, ‘Events’ and ‘CensoringOne’ should be specified appropriately, similarly to fit_tte_data or event_prediction. A series of arguments governing the appearance of the plot are then available; see the individual help file for plot_km_fit for more details.
The plot displays the Kaplan Meier curve, the fitted event curve and the specified confidence interval. If a fit_KM input is provided, no confidence interval will be displayed as it is not available.
The variance function for the fitted trajectory is calculated from the covariance matrix and parameter estimates using the delta method. The interval is then estimated using moment-fitted beta distributions (mean and variance). These intervals are asymmetric and have better coverage and properties than using a normal approximation, particularly where the variance is large or for estimates close to 0 or 1.
In general, for event prediction to be effective, it is necessary for a close fit between the observed data and the fitted curve. This is firstly due to the obvious need to model closely the behaviour of patients for times that have already been observed, but secondly also because fitting the observed data well is necessary in order to have a plausible chance of extrapolating to future times. The key exception to this principle is where overfitting occurs, but this ought not to be relevant for distributions with only one or two parameters.
3.10 Example
We first create a new example data set with piecewise recruitment using the simulate_trials function in gestate:
# Create the full length trial showing what 'actually' happens
trial_long <- simulate_trials(active_ecurve=Weibull(50,0.8),control_ecurve=Weibull(50,0.8),rcurve=recruit,fix_events=200,iterations=1,seed=12345,detailed_output=TRUE)
# Create a shortened version with data for use in the event prediction
trial_short <- set_assess_time(data=trial_long,time=10,detailed_output = FALSE)
# Create the trajectory of events
maxtime <- max(ceiling(trial_long[,"Assess"]))
events <- rep(NA,maxtime)
for (i in 1:maxtime){
events[i] <- sum(1-set_assess_time(trial_long,i)[,"Censored"])
}The ‘trial_short’ object contains the patient data observable at the time of event prediction (10 months). Note that a ‘1’ in the Censored column denotes a censoring, so the censoringOne argument needs to be set to TRUE, and that the unit of time is in months, so the units argument must be “Months”. For real patient data it is more usual to have units of days, which is the default.
‘recruit’ (see section 2.3) contains the observed and predicted patient recruitment as an RCurve object. ‘events’ contains the numbers of events observed at each time point (both before and after prediction) and can be used to assess performance. At the time of event prediction, 10, there are 47 observed events and the number at risk is 500-47= 453; these will be the conditioning parameters. We assume that there is no dropout censoring.
We can now perform a simple event prediction using Weibull curve type selection:
We can then plot the conditional trajectory, with prediction intervals, alongside the events that were (or will be) observed at each time point:
# Plot observed events and conditional predictions
plot_ep(predictions,trajectory="conditional",which_PI="conditional",max_time=40,observed=events,target=200,max_E=200)
The plot shows a reasonable agreement between the observed trajectory and the estimated conditional trajectory, which predicts 200 events would be reached by month 29. However, the 95% prediction interval is wide (particularly on the lower bound), and the prediction interval for the time to reach 200 events is 22-55 months. In practice therefore, event prediction at this time point is not that useful for estimating a trial end date.
We can now revisit the event prediction at two later time points (14 and 18 months), and compare the prediction accuracies:
trial_less_short <- set_assess_time(data=trial_long,time=14,detailed_output = FALSE)
predictions2 <- event_prediction(data=trial_less_short, Event="Censored",censoringOne=TRUE, type="Weibull", rcurve=recruit,max_time=60, cond_Events=101,cond_NatRisk=399,cond_Time=14, units="Months")
trial_mature <- set_assess_time(data=trial_long,time=18,detailed_output = FALSE)
predictions3 <- event_prediction(data=trial_mature, Event="Censored",censoringOne=TRUE, type="Weibull", rcurve=recruit,max_time=60, cond_Events=148,cond_NatRisk=352,cond_Time=18, units="Months")
Only 4 months later, the 95% prediction interval is down to 23-34 months, with an estimate of 27 months, while by 18 months, a tight 95% prediction interval of 23-27 months is possible around an estimate of 25 months. The tightening of the prediction interval over time is partly due to a reduction in the number of events remaining to be predicted, but also the increasing precision with which the parameters of the parametric Kaplan Meier curve fit. As a rule of thumb, 100 observed events appear to be required to give reasonable precision of estimation of Weibull parameters, but the greater the length of extrapolation required, the more precise the estimates have to be to give useful accuracy.
The conditional prediction interval can therefore demonstrate that in this case where a target of 200 events is desired, event prediction is probably too inaccurate to be of use at 10 months, somewhat indicative at 14 months, and fairly accurate by 18 months. As the trial has ‘progressed’, the estimate of the time of finishing has dropped from 29 months to 27 months and then again to 25 months. Depending on tolerance for risk, it may be more useful to use a smaller prediction interval e.g. 80 or 90%.
It should be noted that all of these prediction intervals are dependent upon the underlying distribution being Weibull, and that the same distribution applies to all patients, irrespective of e.g. the time of recruitment. They also depend on the recruitment distribution being fixed, i.e. containing no variability - this may be unreasonable for predictions at a very early time point. These limitations will tend to lead to the interval being narrower than it ought. Conversely, the prediction intervals do not account for knowledge obtained from prior trials or beliefs from the trial design stage; if there are strong prior expectations aligned with the observed data then the prediction intervals may be wider than they need be.
Finally, we can check the quality of the curve fitting and the degree of uncertainty around it.
From visual inspection, in all three cases there is reasonable agreement between the KM and fitted curves (albeit there is not much data at 10 months). This ensures their predictions are likely to at least be plausible. The width of the confidence interval using 10 months data is roughly double that for 18 months data. It is clear that the 10-month time point gives very little confidence about the Survival function when extrapolating. However, the 14 and 18 month time points provide some information about likely future behaviour, so do appear to be potentially useful for event prediction as long as the assumption of a Weibull distribution is reasonable.