General-Purpose Methods for Simulating Survival Data for Expected Value of Sample Information Calculations

Abstract

Background

Expected value of sample information (EVSI) quantifies the expected value to a decision maker of reducing uncertainty by collecting additional data. EVSI calculations require simulating plausible data sets, typically achieved by evaluating quantile functions at random uniform numbers using standard inverse transform sampling (ITS). This is straightforward when closed-form expressions for the quantile function are available, such as for standard parametric survival models, but these are often unavailable when assuming treatment effect waning and for flexible survival models. In these circumstances, the standard ITS method could be implemented by numerically evaluating the quantile functions at each iteration in a probabilistic analysis, but this greatly increases the computational burden. Thus, our study aims to develop general-purpose methods that standardize and reduce the computational burden of the EVSI data-simulation step for survival data.

Methods

We developed a discrete sampling method and an interpolated ITS method for simulating survival data from a probabilistic sample of survival probabilities over discrete time units. We compared the general-purpose and standard ITS methods using an illustrative partitioned survival model with and without adjustment for treatment effect waning.

Results

The discrete sampling and interpolated ITS methods agree closely with the standard ITS method, with the added benefit of a greatly reduced computational cost in the scenario with adjustment for treatment effect waning.

Conclusions

We present general-purpose methods for simulating survival data from a probabilistic sample of survival probabilities that greatly reduce the computational burden of the EVSI data-simulation step when we assume treatment effect waning or use flexible survival models. The implementation of our data-simulation methods is identical across all possible survival models and can easily be automated from standard probabilistic decision analyses.

Highlights

Expected value of sample information (EVSI) quantifies the expected value to a decision maker of reducing uncertainty through a given data collection exercise, such as a randomized clinical trial. In this article, we address the problem of computing EVSI when we assume treatment effect waning or use flexible survival models, by developing general-purpose methods that standardize and reduce the computational burden of the EVSI data-generation step for survival data.

We developed 2 methods for simulating survival data from a probabilistic sample of survival probabilities over discrete time units, a discrete sampling method and an interpolated inverse transform sampling method, which can be combined with a recently proposed nonparametric EVSI method to accurately estimate EVSI for collecting survival data.

Our general-purpose data-simulation methods greatly reduce the computational burden of the EVSI data-simulation step when we assume treatment effect waning or use flexible survival models. The implementation of our data-simulation methods is identical across all possible survival models and can therefore easily be automated from standard probabilistic decision analyses.

Keywords

economic evaluation model expected value of sample information simulation methods survival data value of information

Expected value of sample information (EVSI) quantifies the expected value to a decision maker of reducing uncertainty through a given data collection exercise, such as a randomized clinical trial.¹ Methods for computing the EVSI for collecting survival data (i.e., time-to-event data) when there is uncertainty about the choice of survival model have recently been developed by Vervaart et al.² These methods require, in common with other EVSI methods, the simulation of plausible study data sets that reflect the study design proposed for collecting future data and the time-to-event distribution of individuals included in such a study.³ This is typically achieved by evaluating quantile functions at random uniform numbers using standard inverse transform sampling (ITS). The standard ITS method is straightforward to implement when closed-form expressions for the quantile function are available, such as for standard parametric survival models, but these are often not available when assuming treatment effect waning and for flexible survival models.

Sufficient evidence on time-to-event outcomes, such as overall survival (OS) and time to progression, is crucial for accurately determining the long-term effects of new treatments.⁴ Yet, health technology assessments often have to rely on immature survival data obtained from trials at an early stage, especially for new cancer treatments.⁵ This can partly be explained by the introduction of accelerated licensing schemes for new pharmaceuticals by regulatory bodies such as the European Medicines Agency^6,7 and the US Food and Drug Administration.⁸ Immature survival data require a high degree of extrapolation, which led to the introduction of flexible survival models such as response-based landmark models, mixture cure models, relative survival models, and model-averaging approaches.^9,10 Nevertheless, more complex models do not necessarily result in plausible extrapolations, and therefore, extrapolations are often supplemented with assumptions about disease progression and treatment mechanisms. For example, the National Institute for Health and Care Excellence recommends that waning of treatment effects is considered in technology appraisals,¹¹ for instance, by assuming no more treatment benefit beyond a chosen time point¹² or by assuming that the treatment effect diminishes over the long term.¹³ This is typically implemented in cost-effectiveness models by adjusting the predicted hazards, thereby altering the survival probabilities generated by parametric survival models. This poses a challenge for the EVSI data-simulation step for survival data, as closed-form expressions for the quantile function are often unavailable for custom distributions that incorporate assumptions about treatment effect waning and for flexible survival models. In these circumstances, the standard ITS method for simulating survival data could be implemented by numerically evaluating the quantile distributions at each iteration in a probabilistic analysis, but this can greatly increase the computational burden.

In this article, we address the problem of computing EVSI when we assume treatment effect waning or use flexible survival models, by developing general-purpose methods that standardize and reduce the computational burden of the EVSI data-generation step for survival data. We develop a discrete sampling method and an interpolated ITS method for simulating survival data from a probabilistic sample of survival probabilities over discrete time units. The discrete sampling method samples time cycles using the survival probabilities and sets the event times to the half-cycle times. The interpolated ITS method extends this to continuous time by initially sampling random uniform numbers between 0 and 1 and then interpolating the survival probabilities using cubic splines at the sampled numbers and recording the interpolated cycle times. We demonstrate in an illustrative case study that, when the general-purpose data-simulation methods are combined with a recently proposed nonparametric EVSI method,^2,14 EVSI computations for survival data can be automated from standard probabilistic decision analyses, irrespective of the assumed data-generating process.

The article is structured as follows. In the second section, we first describe the standard ITS method and then introduce the general-purpose methods for simulating survival data. In the third section, we compare the standard ITS method and the general-purpose data-simulation methods based on an illustrative partitioned survival model for scenarios with and without adjustment for treatment effect waning. In the final section, we conclude with a short discussion.

Method

Decision Problem

In health technology assessment, cost-effectiveness models are widely used to compare alternative health technologies in terms of expected costs $V$ and health benefits $Q$ relative to a cost-effectiveness threshold $λ$ .¹⁵ There is usually a choice between a small number of $D$ decision options, indexed $d = 1, \dots D$ , such as a new treatment and standard care. The cost-effectiveness of decision option $d$ can be expressed in terms of net health benefit, ${NB}_{d} = Q_{d} - V_{d} / λ$ , or net monetary benefit, ${NB}_{d} = Q_{d} λ - V_{d}$ . A cost-effectiveness model, which we denote as ${NB}_{d} (θ)$ , predicts the net benefit of decision option $d$ , given a vector of model input parameters $θ$ , such as probabilities, costs, and health-related quality-of-life weights. We represent our current knowledge about the cost-effectiveness model parameters in the joint probability distribution $p (θ)$ . The joint distribution $p (θ)$ is typically defined using a combination of statistical models fitted to individual patient data, external data sources, clinical expert opinion, and assumptions about biological plausibility, and there may also be dependency between elements of $θ$ .

EVSI for a New Study

The expected value of the optimal decision given current information is the value of the decision option that maximizes expected net benefit and can therefore be considered cost-effective,

\max_{d} E_{θ} {{NB}_{d} (θ)},

(1)

A new study would provide data $X$ relating to $θ$ , which could for example be a randomized controlled trial that collects OS data for a new treatment and standard care. Once the study data have been collected, we update the joint distribution of the cost-effectiveness model parameters with the new data using Bayes theorem, giving a posterior distribution $p (θ | X)$ .

The expected value of the optimal decision made after collecting $X$ is

\max_{d} E_{θ | X} {{NB}_{d} (θ)} .

(2)

The value of collecting $X$ is derived from its potential to help a decision maker avoid recommending a treatment that reduces net benefit, as we may learn that the decision option that is considered cost-effective given current information turns out to be cost-ineffective given new information. However, as we have not yet collected $X$ , we must take the expectation over the distribution of all possible data sets, $p (X)$ . This requires simulating plausible data sets from the distribution of the new data, $x ~ p (X)$ . To simulate $x$ , we usually need to specify a parametric data-generating distribution, $p (X | θ)$ , which depends on the type of data that will be generated by the proposed study and how it will be used to update the model parameters.³ We can simulate data sets from $p (X)$ by first sampling from the joint distribution of the model parameters, $θ^{*} ~ p (θ)$ , and then sampling from the parametric data-generating distribution given the parameter sample, $x^{*} ~ p (X | θ^{*})$ . This gives us a pair of samples ${x^{*}, θ^{*}}$ from the joint distribution $p (X, θ)$ , and therefore the samples $x^{*}$ are drawn from the marginal distribution $p (X)$ .

The expected value of the decision made with additional sample information is given by

E_{X} [\max_{d} E_{θ | X} {{NB}_{d} (θ)}] .

(3)

The EVSI¹⁶ measures the expected value of reducing uncertainty about the optimal decision by collecting $X$ , which is defined as the difference between equation (3) and equation (1),

EVSI = E_{X} [\max_{d} E_{θ | X} {{NB}_{d} (θ)}] - \max_{d} E_{θ} {{NB}_{d} (θ)} .

(4)

In the next section, we describe the EVSI estimation procedure for collecting time-to-event data.

Computing EVSI for Time-to-Event Data

Time-to-Event Data

Time-to-event data, such as time to disease progression and time to death, are frequently collected in the context of clinical trials. A special feature of time-to-event data is censoring, which occurs when the follow-up time is not long enough to observe the event of interest for all individuals or when individuals are lost to follow-up.¹⁷ A single time-to-event data set $x$ , collected between time 0 and $t_{o b s}$ , consists of $i = 1, \dots, n$ survival times $x_{i}$ and censoring indicators $δ_{i}$ , $x = {x_{1}, \dots, x_{n}, δ_{1}, \dots, δ_{n}}$ , where $δ_{i} = 1$ when $x_{i}$ is an observed event and $δ_{i} = 0$ when $x_{i}$ is a censored observation.

To predict outcomes over the long term, censored time-to-event data usually need to be extrapolated beyond the observed follow-up period using a parametric survival model.⁴ Parametric models are commonly specified using either the survivor function, $S (t, θ)$ , or hazard function, $h (t, θ)$ , with time denoted as $t$ .

The survivor function $S (t)$ defines the probability of survival up to time $t$ , given by

S (t, θ) = P (T > t) = 1 - F (t, θ), 0 < t < inf,

(5)

where $F (t, θ)$ is the cumulative distribution function.

The hazard function $h (t, θ)$ defines the instantaneous event rate at time $t$ conditional on survival up to time $t$ ,

h (t, θ) = lim_{δ t \to 0} {\frac{P (t \leq T < t + δ t | T \geq t)}{δ t}} = \frac{f (t, θ)}{S (t, θ)},

(6)

where $f (t, θ)$ is the probability density function.

Most cost-effectiveness models are in discrete time and therefore evaluate $S (t, θ)$ at fixed time intervals based on a set amount of time called model cycles, thereby generating a vector of survival probabilities over discrete time units for each sampled value for $θ$ , from which the distribution of patients across the health states over time can be constructed. The expected costs, life-years, and quality-adjusted life-years of a new technology can be estimated by integrating the distribution of patients across the health states and the health state values, which reflect the costs and health-related quality-of-life weights associated with each health state.¹⁵

Computing the Expected Net Benefits Given Current Information

We can compute the expected net benefits given current information in a probabilistic analysis (PA) using Monte Carlo simulation. This involves sampling $k = 1, \dots, K$ values, $θ^{(k)}$ , from the distribution of the model parameters, $p (θ)$ , and then evaluating the survivor functions, $S (t_{c}, θ^{(k)})$ , given model cycle $c = 0, \dots, C$ . This results in $K$ vectors of survival probabilities, $s^{(k)} = {s_{0}^{(k)}, \dots, s_{C}^{(k)}}$ , where $s_{c}^{(k)} = S (t_{c}, θ^{(k)})$ . We assume $s \subset θ$ and compute the net benefits conditional on the survival curves and other model parameters, ${NB}_{d} (θ^{(k)})$ , for each $d$ . We compute the expected net benefits given current information by averaging over the $K$ net benefits for each $d$ .

Table 1 illustrates a PA sample in which $K$ vectors of survival probabilities $s^{(k)}$ have been generated for a single decision option $d$ .

Table 1

Probabilistic Analysis Sample with $K$ Sampled Vectors of Survival Probabilities $s^{(k)}$ for a Single Decision Option

Model Cycle	$s^{(1)}$	$s^{(2)}$	$s^{(3)}$	$\dots$	$s^{(K)}$
0	1.00	1.00	1.00	$\dots$	1.00
1	0.85	0.79	0.88	$\dots$	0.92
2	0.76	0.70	0.79	$\dots$	0.88
3	0.69	0.63	0.72	$\dots$	0.79
$⋮$	$⋮$	$⋮$	$⋮$	$⋱$	$⋮$
$C$	0.00	0.00	0.00	$\dots$	0.00
$\sum_{c = 0}^{C} s_{c}^{(k)}$	14.49	10.82	12.11	$\dots$	18.79

Simulating Time-to-Event Data

Standard ITS Method

To compute the expected net benefits given new time-to-event data $x$ , we must first sample plausible data sets from the distribution of the new data conditional on the sampled parameter values, $x^{(k)} ~ p (X | θ^{(k)})$ . We can simulate random survival times from a parametric survival distribution using standard functions, such as the rweibull function in R.¹⁸ These standard functions use ITS, which is a method for generating random values from an arbitrary distribution. Let $F$ be a continuous cumulative distribution function (CDF), with inverse CDF (i.e., quantile function) denoted by $F^{- 1}$ . The intuition behind the ITS method is that if we define $T = F^{- 1} (U)$ , where $U ~ Unif [0, 1]$ , then $T$ is distributed according to $F$ , that is, $F (T)$ . Since the survivor function is defined as $S (t, θ) = 1 - F (t, θ)$ , we can generate survival times by sampling values $u_{i}, i = 1, \dots, n$ , from a uniform distribution on the interval $[0, 1]$ and plugging these into the inverse survivor function, $S^{- 1} (u_{i}, θ)$ . The inverse survivor function can be derived from the inverse CDF, $S^{- 1} (1 - u_{i}, θ)$ = $F^{- 1} (u_{i}, θ)$ , and if $U ~ Unif [0, 1]$ , it follows that $1 - U$ is also $~ Unif [0, 1]$ . We then censor the survival times at the new follow-up time $t_{n e w}$ .

The standard ITS scheme for simulating time-to-event data is given in Box 1.

Box 1

Standard Inverse Transform Sampling Scheme for Simulating Time-to-Event Data

for

k = 1, \dots, K

parameter samples

θ^{(k)}

do
Sample a vector of values

u^{(k)} = {u_{1}^{(k)}, \dots, u_{n}^{(k)}}

, with each element drawn from

Unif [0, 1]

Generate a dataset of

n

survival times

x^{(k)} = {x_{1}^{(k)}, \dots, x_{n}^{(k)}}

by plugging

u^{(k)}

into

S^{- 1} (u^{(k)}, θ^{(k)})

Censor each element of

x^{(k)}

at the new follow-up time

t_{n e w}

end

Sampling from a Weibull distribution, such as by using the rweibull function in R, can result in any value between 0 and infinity. If we want to ensure that the sampled survival times do not exceed a biologically plausible time horizon $t_{h}$ , which would typically be equal to the time horizon used in the decision model, then we need to sample $u_{i}^{(k)}$ on the interval $[p_{1}^{(k)}, 1]$ , where $p_{1}^{(k)}$ is the survivor function evaluated at $t_{h}$ , that is, $S (t_{h}, θ^{(k)})$ . If, however, the chance of simulating implausible survival times that exceed $t_{h}$ is substantial, the analyst should consider using a more realistic survival model, in particular models that incorporate background mortality and other external information.^9,19 However, this evidence may not be available and is exactly the evidence that is under consideration to collect in the value-of-information analysis.

The algorithm in Box 1 could, in theory, be used for any survival model for which we can define a hazard function, $h (t, θ)$ . The hazard could be integrated to produce a cumulative hazard function, $H (t, θ)$ , from which we can derive the survivor function, $S (t, θ) = \exp [- H (t, θ)]$ . The survivor function could then be inverted to produce the quantile function, $S^{- 1} (\cdot, θ)$ , and survival times could be simulated by evaluating $S^{- 1} (\cdot, θ)$ at a random uniform sample, as described in Box 1.

Analytic solutions to the integrals and function inverses, such as implemented in the rweibull function, may not be available for flexible survival models, such as relative survival models, spline models, mixture cure models, and response-based landmark models^9,10 and for custom distributions, including hazard functions that incorporate assumptions about treatment effect waning. In these circumstances, the integrals and function inverses could be evaluated numerically, for example, by using the integrate() and uniroot() functions in base R. The flexsurv and msm packages on CRAN also have a function qgeneric() designed to invert a generic CDF. These numerical methods, however, greatly increase the computational burden of the EVSI computations, as the integration and inversion steps need to be repeated for each of $k = 1, \dots, K$ simulations. In Appendix A, we provide a simple example and step-by-step implementation in R of the standard ITS method based on analytic and numerical solutions. The example illustrates that even for a simple exponential model, numerically evaluating the integrals and function inverses results in a greater than 10,000-fold increase in computational time (Table A1).

In the next section, we will introduce general-purpose methods for simulating survival data that can be standardized from standard probabilistic decision analyses and greatly reduce the computational burden of the standard ITS method when closed-form expressions for the quantile function are unavailable.

Simulating Time-to-Event Data from a Vector of Survival Probabilities over Discrete Time Units

Interpolated ITS method

We can also use the ITS method to generate $x$ from a vector of survival probabilities over discrete time units, as illustrated in Figure 1. Our PA sample consists of $K$ vectors of survival probabilities, $s^{(k)} = {s_{0}^{(k)}, \dots, s_{C}^{(k)}}$ , given model cycle $c = 0, \dots, C$ . We can approximate values from the quantile function by interpolating the vectors of survival probabilities, $S^{- 1} (\cdot, s^{(k)})$ . We can achieve this by fitting monotone cubic splines, a type of piecewise polynomial interpolation, to each consecutive set of survival probabilities and cycle times, ${(s_{0}^{(k)}, s_{1}^{(k)}), (t_{0}, t_{1})}, \dots, {(s_{C - 1}^{(k)}, s_{C}^{(k)}), (t_{C - 1}, t_{C})}$ . We can then sample survival times from the approximate quantile function $S^{- 1} (\cdot, s^{(k)})$ using the ITS method.

Figure 1

Illustration of the interpolated inverse transform sampling method for simulating time-to-event data. A survival time of 3.6 y has been simulated by first sampling a value of 0.32 from a uniform distribution between 0 and 1 and then interpolating the survival probabilities over discrete time units using monotone cubic splines at 0.32 and recording the interpolated cycle time of 3.6 y.

The interpolated ITS sampling scheme for simulating time-to-event data from a vector of survival probabilities over discrete time units is given in Box 2.

Box 2

Interpolated Inverse Transform Sampling Scheme for Simulating Time-to-Event Data from a Vector of Survival Probabilities over Discrete Time Units

for

k = 1, \dots, K

vectors of survival probabilities

s^{(k)} = {s_{0}^{(k)}, \dots, s_{C}^{(k)}}

do
Sample a vector of values

u^{(k)} = {u_{1}^{(k)}, \dots, u_{n}^{(k)}}

, with each element drawn from

Unif [s_{C}^{(k)}, 1]

Generate a dataset of

n

survival times

x^{(k)} = {x_{1}^{(k)}, \dots, x_{n}^{(k)}}

by interpolating

s^{(k)}

u^{(k)}

using monotone cubic splines and recording the interpolated cycle times
Censor each element of

x^{(k)}

at the new follow-up time

t_{n e w}

end

Discrete sampling method

An alternative approach for simulating $x$ from a vector of survival probabilities over discrete time units $s^{(k)}$ is by sampling discrete “bins” of cycle times with probability equal to the cumulative density of each bin, which is illustrated in Figure 2. We set the value of each bin equal to the half-cycle time $t_{c + 0.5}$ , $c = 0, \dots, C - 1$ , that is, the midpoint in between each consecutive set of cycle times. We can estimate the cumulative density of each bin from the survival probabilities, $(s_{c}^{(k)} - s_{c + 1}^{(k)})$ . We can then generate survival times by sampling ${t_{0.5}, \dots, t_{C - 0.5}}$ with probability ${(s_{0}^{(k)} - s_{1}^{(k)}), \dots, (s_{C - 1}^{(k)} - s_{C}^{(k)})}$ .

Figure 2

Illustration of the discrete sampling method for simulating time-to-event data. Survival times can be simulated by sampling from the half-cycle times on the x-axis with probability derived from the survival probabilities over discrete time units, as indicated by the arrows.

The discrete sampling scheme for simulating time-to-event data from a vector of survival probabilities over discrete time units is given in Box 3.

Box 3

Discrete sampling scheme for simulating time-to-event data from a vector of survival probabilities over discrete time units.

for

k = 1, \dots, K

vectors of survival probabilities

s^{(k)} = {s_{0}^{(k)}, \dots, s_{C}^{(k)}}

do
Generate a dataset of

n

survival times

x^{(k)} = {x_{1}^{(k)}, \dots, x_{n}^{(k)}}

by sampling

n

values from the half-cycle times,

{t_{0.5}, \dots, t_{C - 0.5}}

, with probability

{(s_{0}^{(k)} - s_{1}^{(k)}), \dots, (s_{C - 1}^{(k)} - s_{C}^{(k)})}

Censor each element of

x^{(k)}

at the new follow-up time

t_{n e w}

end

Appendix A describes a step-by-step implementation in R of the interpolated ITS and discrete sampling methods, and a comparison of their computational efficiency with the standard ITS method based on analytic and numerical solutions.

Computing the expected net benefits given new data

We could compute the expected net benefits given the $k = 1, \dots, K$ sampled data sets for each decision option, $x_{d}^{(k)}$ , using a nested Monte Carlo scheme, but this can be very computationally demanding as it requires sampling a large number of values from the posterior distribution of the model parameters conditional on each simulated data set, $p (θ | x_{d}^{(k)})$ . A number of efficient approximation methods have therefore been developed in recent years that reduce the computational burden of the nested Monte Carlo procedure to EVSI.²⁰ One of these efficient EVSI methods is a nonparametric regression-based method developed by Strong et al.¹⁴ that does not require a parametric distribution for the data. The regression-based method for computing EVSI relies on estimating the functional relationship between the posterior expected net benefits and the simulated data sets, thereby avoiding the need to sample from the posterior distributions $p (θ | x_{d}^{(k)})$ , as is required for the nested Monte Carlo scheme.

In the regression-based approach, we require only the vectors of $k = 1, \dots, K$ prior net benefits ${NB}_{d} (θ^{(k)})$ for each decision option $d$ that we generated in the PA, and the corresponding data sets $x_{d}^{(k)}$ that we have simulated using either of the sampling schemes above. The observed net benefits ${NB}_{d} (θ^{(k)})$ can be expressed as a sum of the conditional expectation of the net benefit given the data, $E_{θ | x_{d}^{(k)}} {{NB}_{d} (θ)}$ , which we require to estimate the EVSI (equation [4]), and a mean-zero error term, $ε^{(k)}$ ,

{NB}_{d} (θ^{(k)}) = E_{θ | x_{d}^{(k)}} {{NB}_{d} (θ)} + ε^{(k)} .

(7)

Strong et al.¹⁴ explain that the conditional expectation $E_{θ | x_{d}^{(k)}} {{NB}_{d} (θ)}$ can be thought of as an unknown function of $x_{d}^{(k)}$ . We denote this function $g (x_{d}^{(k)})$ and substitute this into equation (7), giving

{NB}_{d} (θ^{(k)}) = g (x_{d}^{(k)}) + ε^{(k)} .

(8)

We then summarize $x_{d}^{(k)}$ using a low-dimensional summary statistic for each $d$ , $T (x_{d}^{(k)})$ ,

{NB}_{d} (θ^{(k)}) = g {T (x_{d}^{(k)})} + ε^{(k)} .

(9)

A convenient choice for $T (x_{d}^{(k)})$ is the number of observed events $e_{d}^{(k)}$ and the total time at risk $y_{d}^{(k)}$ for each simulated data set $x_{d}^{(k)}$ , that is, $T (x_{d}^{(k)}) = {e_{d}^{(k)}, y_{d}^{(k)}}$ , which has been shown to give good results for various survival models.²

We can estimate the posterior net benefits by regressing the prior net benefits, ${NB}_{d} (θ^{(k)})$ , on the summary statistic $T (x_{d}^{(k)})$ ,

{NB}_{d} (θ^{(k)}) = g {T (x_{d}^{(k)})} + ε,

(10)

where $g_{d}$ is a function of the data $x_{d}^{(k)}$ for each $d$ and $ε$ is an error term with zero mean. We can achieve this by fitting a generalized additive model (GAM), which is a flexible nonparametric regression model, to each $d$ and extracting the regression model fitted values ${\hat{g}}_{d}^{(k)}$ , which are estimates of the posterior net benefits.

The GAM-based EVSI estimate is given by

EVSI ≃ \frac{1}{K} \sum_{k = 1}^{K} \max_{d} {\hat{g}}_{d}^{(k)} - \max_{d} \frac{1}{K} \sum_{k = 1}^{K} {\hat{g}}_{d}^{(k)} .

(11)

EVSI for an ongoing study

When a trial is ongoing at the point of decision making, there could be value in reducing uncertainty by collecting additional data from the ongoing trial before making an adoption decision. This is especially common for cost-effectiveness analyses of new cancer drugs, which increasingly rely on immature data obtained from trials in an early stage.⁵

We denote the new data collected between the observed follow-up time $t_{o b s}$ and future point $t_{n e w}$ as $\tilde{x} = {{\tilde{x}}_{1}, \dots, {\tilde{x}}_{\tilde{n}}, {\tilde{δ}}_{1}, \dots, {\tilde{δ}}_{\tilde{n}}}$ , where $\tilde{n}$ is the number of study participants at risk at $t_{o b s}$ . As we have not yet collected the new data, we need to simulate plausible new data sets conditional on the current data, $\tilde{x} ~ (\tilde{X} | X)$ . The value of extending an existing trial’s follow-up from $t_{o b s}$ to future point $t_{n e w}$ is given by²

\begin{matrix} EVSI (ongoing study) & = E_{\tilde{X} | X} [\max_{d} E_{θ | X, \tilde{X}} {{NB}_{d} (θ)}] \\ - \max_{d} E_{θ | X} {{NB}_{d} (θ)}, \end{matrix}

(12)

where the first term is the expected value of a decision based on the joint posterior distribution of $θ$ given both new data, $\tilde{X}$ , collected between $t_{o b s}$ and $t_{n e w}$ , and current data, $X$ , collected between time zero and $t_{o b s}$ , which is computed by averaging over the posterior net benefits of the decision option that maximizes this quantity conditional on both new data and current data. The second term is the expected value of a decision based on the joint distribution of $θ$ given current data collected up until $t_{o b s}$ .

Events beyond $t_{o b s}$ are conditional on survival up to $t_{o b s}$ . Therefore, to simulate future survival times for $\tilde{n}$ patients at risk at $t_{o b s}$ , we need to sample from a conditional distribution that is left truncated at $t_{o b s}$ . For the standard ITS scheme in Box 1, this requires evaluating the survivor function at $t_{o b s}, p_{2}^{(k)} = S (t_{o b s}, θ^{(k)})$ and then sampling random uniform numbers, $u^{(k)}$ , on the interval $[p_{1}^{(k)}, p_{2}^{(k)}]$ , where $p_{1}^{(k)}$ is the survivor function evaluated at the model time horizon, $S (t_{h}, θ^{(k)})$ . For the interpolated ITS scheme in Box 2, we can find $p_{2}^{(k)}$ by interpolating $s^{(k)}$ at $t_{o b s}$ and then sampling $u^{(k)}$ on the interval $[s_{C}^{(k)}, p_{2}^{(k)}]$ . For the discrete sampling scheme in Box 4, we need to sample from the subset of half-cycle times that are larger than $t_{o b s}$ , i.e. ${t_{c + 0.5}, \dots, t_{C - 0.5}}$ for $t_{c + 0.5} > t_{o b s}$ .

Model-averaged EVSI

Uncertainty about the choice of survival model is often a key driver of decision uncertainty, particularly when data are immature.²¹ If we are uncertain about choosing from a set of competing survival models for extrapolating study data over the long term, $M = M_{r}, r = 1, \dots, R$ , then we could account for this in the EVSI computations² by using model averaging.^22–24 Before we collect new data $X$ , our beliefs about the plausibility of each model is represented by the prior model probabilities $P (M_{r})$ . These could, for example, be derived from the Akaike’s information criterion²⁵ or other measures of model fit and parsimony. In the PA, we then sample a survival model $M_{r}^{(k)}$ with probability $P (M_{r})$ before sampling $θ_{r}^{(k)}$ from the distribution of the parameters of the sampled survival model, $p (θ_{r}, M_{r}^{(k)})$ . Since this changes only the values for $s^{(k)}$ , the discrete sampling scheme in Box 2 and interpolated ITS scheme in Box 3, as well as the GAM-based EVSI estimation following equation (11), are identical to the single-model case.

Synthetic Case Study

Decision Problem and Model Definition

To demonstrate the application of our methods, we developed a simple yet realistic synthetic case study based on a partitioned survival model (PSM)²⁶ comparing a new treatment (d = 1) with standard care (d = 2). The PSM uses OS and progression-free survival (PFS) curves to estimate the proportion of patients in 3 health states: PFS, postprogression survival (PPS), and death, given $c = 0, \dots, 360$ monthly model cycles corresponding to an overall time horizon $t_{h} = 30$ y. We assumed OS and PFS follow independent Weibull distributions for each $d$ , parameterized in terms of log shape $α$ and log scale $β$ . The Weibull model parameters for OS are $θ_{o s 1} = (α_{o s 1}, β_{o s 1})$ for the new treatment and $θ_{o s 2} = (α_{o s 2}, β_{o s 2})$ for standard care, and $θ_{p f s 1} = (α_{p f s 1}, β_{p f s 1})$ and $θ_{p f s 2} = (α_{p f s 2}, β_{p f s 2})$ for PFS, respectively. We estimated the Weibull model parameters using maximum likelihood from a synthetic data set containing 100 OS times and 100 PFS times for each trial arm with a maximum follow-up of 24 mo. Further details about the synthetic case study data set are given in Appendix B. The other model parameters are utility for PFS ( $U_{p f s}$ ), utility for postprogression ( $U_{pps}$ ), drug costs for the new treatment ( $C_{drug 1}$ ) in PFS, medical costs for the new treatment ( $C_{med 1}$ ), medical costs for standard care ( $C_{med 2}$ ), annual discount rate ( $r$ ), and monetary value of 1 quality-adjusted life-year, ( $λ$ ). The definitions and prior distributions for the case study model parameters are given in Table 2, and the net benefit functions are given in Appendix C.

Table 2

Prior Parameter Distributions for the Partitioned Survival Model Parameters

Weibull Survival Model Parameters		$Mean$ , $μ$	Covariance Matrix, $Σ$	Distribution
Overall survival
Log shape for new treatment	$(\begin{matrix} α_{o s 1} \\ β_{o s 1} \end{matrix})$	$(\begin{matrix} 0.312 \\ 4.089 \end{matrix})$	$(\begin{matrix} 0.037 & - 0.036 \\ - 0.036 & 0.056 \end{matrix})$	Bivariate normal $(μ, Σ)$
Log scale for new treatment	$(\begin{matrix} α_{o s 1} \\ β_{o s 1} \end{matrix})$	$(\begin{matrix} 0.312 \\ 4.089 \end{matrix})$		Bivariate normal $(μ, Σ)$
Log shape for standard care	$(\begin{matrix} α_{o s 2} \\ β_{o s 2} \end{matrix})$	$(\begin{matrix} 0.361 \\ 3.842 \end{matrix})$	$(\begin{matrix} 0.029 & - 0.021 \\ - 0.021 & 0.030 \end{matrix})$	Bivariate normal $(μ, Σ)$
Log scale for standard care	$(\begin{matrix} α_{o s 2} \\ β_{o s 2} \end{matrix})$	$(\begin{matrix} 0.361 \\ 3.842 \end{matrix})$		Bivariate normal $(μ, Σ)$
Progression-free survival
Log shape for new treatment	$(\begin{matrix} α_{p f s 1} \\ β_{p f s 1} \end{matrix})$	$(\begin{matrix} 0.161 \\ 3.590 \end{matrix})$	$(\begin{matrix} 0.019 & - 0.010 \\ - 0.010 & 0.021 \end{matrix})$	Bivariate normal $(μ, Σ)$
Log scale for new treatment	$(\begin{matrix} α_{p f s 1} \\ β_{p f s 1} \end{matrix})$	$(\begin{matrix} 0.161 \\ 3.590 \end{matrix})$		Bivariate normal $(μ, Σ)$
Log shape for standard care	$(\begin{matrix} α_{p f s 2} \\ β_{p f s 2} \end{matrix})$	$(\begin{matrix} 0.209 \\ 3.294 \end{matrix})$	$(\begin{matrix} 0.014 & 0.004 \\ - 0.004 & 0.013 \end{matrix})$	Bivariate normal $(μ, Σ)$
Log scale for standard care	$(\begin{matrix} α_{p f s 2} \\ β_{p f s 2} \end{matrix})$	$(\begin{matrix} 0.209 \\ 3.294 \end{matrix})$		Bivariate normal $(μ, Σ)$
Utility parameters		Mean, $μ$	Standard error, SE	Distribution
Utility for progression-free survival	$U_{p f s}$	0.80	0.04	Beta(80,20)
Utility for post-progression survival	$U_{pps}$	0.50	0.05	Beta(50,50)
Cost parameters (monthly)		Mean $μ$	Standard error, SE	Distribution
Drug costs for new treatment	$C_{drug 1}$	1, 200	—	Constant
Medical costs for new treatment	$C_{med 1}$	500	250	Gamma(4,0.008)
Medical costs for standard care	$C_{med 2}$	500	250	Gamma(4,0.008)
Other parameters		Mean, $μ$	Standard error, SE	Distribution
Annual discount rate	$r$	0.035	—	Constant
Monetary value of 1 quality-adjusted life-year	$λ$	80, 000	—	Constant

Treatment-stopping rule and treatment effect waning

We also considered a scenario with a 2-y treatment-stopping rule, after which the drug costs for the new treatment $C_{drug 1} = 0$ . We assumed that, on expectation, the treatment effect on OS and PFS would wane after stopping treatment at year 2 until there was no more treatment effect by year 4. A treatment effect waning is commonly implemented in cost-effectiveness models by setting the hazards in the treatment arm and comparator arm equal at a chosen time point. This approach has a number of limitations. First, this could lead to counterintuitive results when the hazard in the comparator arm is below the hazard in the treatment arm after the treatment duration cutoff, in which case the treatment effect increases during the waning period.¹² Second, setting the hazards in the treatment and comparator arm equal underestimates uncertainty about independent survival endpoints. Third, this approach ignores uncertainty about the start time and duration of the waning period.

To avoid these limitations, we used an alternative approach to implement treatment effect waning. For $k = 1, \dots, K$ , we first sampled values, $t_{w 1}^{(k)}$ and $t_{w 2}^{(k)}$ , for the start time and duration of the waning period from $p (t_{w 1}) ~ LogNormal (2.83, 0.69)$ and $p (t_{w 2}) ~ LogNormal (2.83, 0.69)$ , corresponding to mean times and standard deviations of 24 months for both $t_{w 1}$ and $t_{w 2}$ . In practical applications, $p (t_{w 1})$ and $p (t_{w 2})$ could be estimated from previous studies or by expert elicitation.²⁷ We then computed vectors of “waning hazards” from the Weibull mean survival probabilities for the new treatment and standard care that in linear proportions increase from no waning to full waning over the time periods defined by $t_{w 1}^{(k)}$ and $t_{w 2}^{(k)}$ , and added these to the sampled Weibull hazards for the new treatment, $h_{1}^{(k)}$ .

The waning-adjusted hazard function for OS and PFS for the new treatment is given by

\begin{matrix} Waning - adjusted hazard h_{1}^{*} (t) = \frac{e^{θ_{α 1}}}{e^{θ_{β 1}}} {(\frac{t}{e^{θ_{β 1}}})}^{e^{θ_{α 1}} - 1} \\ + F_{U n i f} (\frac{t - t_{w 1}}{t_{w 2}} | 0, 1) \max [0, f_{s p l i n e} (t, {\hat{h}}_{2} - {\hat{h}}_{1})], \end{matrix}

(13)

where the first term is the Weibull hazard function, $F_{U n i f} (\cdot | 0, 1)$ is a uniform cumulative distribution function on the interval $[0, 1]$ evaluated at $(t - t_{w 1}) / t_{w 2}$ , and $f_{s p l i n e} (t, {\hat{h}}_{2} - {\hat{h}}_{1})$ is a spline function fitted to the difference in mean survival hazards for the new treatment and standard care, ${\hat{h}}_{1}$ and ${\hat{h}}_{2}$ , respectively, with mean survival hazards for treatment option $d$ given by

Mean survival hazards {\hat{h}}_{d} = - \log (1 - \frac{{{\bar{s}}_{d, 1}, \dots, {\bar{s}}_{d, C}}}{{{\bar{s}}_{d, 0}, \dots, {\bar{s}}_{d, C - 1}}}),

(14)

where ${\bar{s}}_{d, c}$ , $c = 0, \dots, C$ , are mean Weibull survival probabilities over discrete time units for each $d$ , given by

Mean survival probability {\bar{s}}_{d, c} = \frac{1}{K} \sum_{k = 1}^{K} e^{- {(t_{c} / e^{β_{d}^{(k)}})}^{e^{α_{d}^{(k)}}}} .

(15)

All other model assumptions and net benefit functions are as above.

The expected Weibull survival curves for the scenarios with and without adjustment for treatment effect waning are given in Figure 3.

Figure 3

Expected Weibull survival curves for overall survival and progression-free survival for the new treatment and for standard care without adjustment (left) and with adjustment for treatment effect waning (right).

Computations

We assume we want to compute the EVSI for a new study that will collect OS and PFS data for the new treatment and for standard care. We considered a sample size of $n = 100$ study participants for each treatment arm and a study follow-up period $t_{n e w} = 1, 2, 3, 4, 5$ and $10$ y. We first performed a PA in which we sampled $K = 2, 000$ parameter values $θ^{(k)}$ from the distribution of the model parameters $p (θ)$ and then evaluated the PSM to obtain the net benefits for each $d$ , ${NB}_{d} (θ^{(k)})$ .

Simulating OS and PFS data

In the scenario without adjustment for treatment effect waning, we simulated OS data sets, $x_{o s 1}^{(k)}$ for the new treatment and $x_{o s 2}^{(k)}$ for standard care, with $n = 100$ times for each treatment arm using the standard ITS method, the interpolated ITS method and the discrete sampling method following the sampling schemes in Box 1, Box 2, and Box 3, respectively. In the scenario with adjustment for treatment effect waning, we implemented the standard ITS method by numerically integrating the waning hazard function and inverting the waning survivor function, as analytic solutions are unavailable. The implementation of the interpolated ITS method and the discrete sampling method is identical to the scenario without treatment effect waning. We summarized each simulated data set using the number of observed OS events, $e_{o s 1}^{(k)}$ and $e_{o s 2}^{(k)}$ , and the total time at risk for OS, $y_{o s 1}^{(k)}$ and $y_{o s 2}^{(k)}$ , that is, $T (x_{o s 1}^{(k)}) = {e_{o s 1}^{(k)}, y_{o s 1}^{(k)}}$ and $T (x_{o s 2}^{(k)}) = {e_{o s 2}^{(k)}, y_{o s 2}^{(k)}}$ .

We simulated PFS data sets for the new treatment and standard care ( $x_{p f s 1}^{(k)}$ and $x_{p f s 2}^{(k)}$ ) using the same sampling schemes as for OS for the scenarios with and without adjustment for treatment effect waning. To prevent double counting of OS data, we censored the simulated PFS times at the proposed follow-up time or at the time point at which, if at all, the PFS curve crosses and is set equal to the OS curve, whichever is soonest. We then computed the number of observed PFS events, $e_{p f s 1}^{(k)}$ and $e_{p f s 2}^{(k)}$ , and the total time at risk for PFS, $y_{p f s 1}^{(k)}$ and $y_{p f s 2}^{(k)}$ , for each simulated PFS data set, that is, $T (x_{p f s 1}^{(k)}) = {e_{p f s 1}^{(k)}, y_{p f s 1}^{(k)}}$ and $T (x_{p f s 2}^{(k)}) = {e_{p f s 2}^{(k)}, y_{p f s 2}^{(k)}}$ .

Computing EVSI via GAM regression

To reduce the number of regression equations¹⁴ and improve the stability of the EVSI computations,²⁸ we used the incremental net benefit (INB), defined as $INB (θ^{(k)}) = {NB}_{1} (θ^{(k)}) - {NB}_{2} (θ^{(k)})$ . We estimated posterior INB by fitting a single GAM model with $INB (θ^{(k)})$ as the dependent variable and separate sets of summary statistics per treatment arm, ${e_{o s 1}^{(k)}, y_{o s 1}^{(k)}, e_{p f s 1}^{(k)}, y_{p f s 1}^{(k)}}$ and ${e_{o s 2}^{(k)}, y_{o s 2}^{(k)}, e_{p f s 2}^{(k)}, y_{p f s 2}^{(k)}}$ , as the independent variables. We implemented the GAM regression using the R package mgcv²⁹ and specified a tensor product cubic regression spline basis for the independent variables, with a maximum basis dimension of 4 to prevent the model from estimating too many coefficients. This has syntax gam (inb∼te(e_os1, y_os1, e_pfs1, y_pfs1, k=4) + te(e_os2, y_os2, e_pfs2, y_pfs2, k=4)). We extracted the GAM model fitted values, ${\hat{g}}^{(k)}$ , which are estimates of the posterior INB, and estimated EVSI using the equation given by

EVSI \approx \frac{1}{K} \sum_{k = 1}^{K} \max {0, {\hat{g}}^{(k)}} - \max \frac{1}{K} \sum_{k = 1}^{K} {0, {\hat{g}}^{(k)}} .

(16)

We computed 95% intervals for the GAM estimator by sampling 2,000 values from a multivariate normal distribution of the GAM coefficients, as described in an appendix of the article by Strong et al.³⁰

Results

Figure 4 shows the EVSI values and 95% intervals without adjustment for treatment effect waning for follow-up times of 1, 2, 3, 4, 5, and 10 y. There is excellent agreement between the standard ITS method, interpolated ITS method, and discrete sampling method both for OS only and OS and PFS. The EVSI reflects diminishing marginal returns for increasing follow-up durations, ranging from 2,711 to 11,311 for OS only and from 3,387 to 11,547 for OS and PFS, and converges toward the partial EVPI for the respective sets of model parameters. This indicates that the value of reducing uncertainty about PFS in addition to OS is relatively small. We did not compute the EVSI for PFS only, since PFS is a composite endpoint that is defined as time to progression or time to death, whichever is soonest, and therefore also requires the collection of OS data. The total computation times for the data-simulation procedures in the scenario without adjustment for treatment effect waning are 12 s, 32 s, and 12 s for the standard ITS method, interpolated ITS method, and discrete sampling method, respectively.

Figure 4

EVSI values for the synthetic case study without adjustment for treatment effect waning. Total computation times for the data-simulation procedures are 12 s (standard ITS), 32 s (interpolated ITS) and 12 s (discrete sampling). EVSI, expected value of sample information; ITS, inverse transform sampling.

The EVSI estimates in the scenario with adjustment for treatment effect waning (Figure 5) are greater than in the scenario without adjustment for treatment effect waning, reflecting the added value of learning about treatment effect waning. The EVSI estimates range from 5,433 to 11,961 for OS only and from 6,064 to 12,384 for OS and PFS, almost twice as high for the 1-y follow-up period compared with the scenario without adjustment for treatment effect waning. The interpolated ITS method and discrete sampling method again agree closely with the standard ITS method but at a greatly reduce computational cost. The computation times for the interpolated ITS and discrete sampling methods are, in fact, the same as in the scenario without adjustment for treatment effect waning, and approximately 3,600 and 10,000 times faster, respectively, than the standard ITS scheme that used numerical solutions for the integrals and function inverses.

Figure 5

EVSI values for the synthetic case study with adjustment for treatment effect waning. Total computation times for the data-simulation procedures are 120,180 s (standard ITS), 33 s (interpolated ITS), and 12 s (discrete sampling). EVSI, expected value of sample information; ITS, inverse transform sampling.

Discussion

Strengths and Limitations

We developed an interpolated ITS method and a discrete sampling method for simulating survival data from a probabilistic sample of survival probabilities over discrete time units. Our general-purpose methods greatly reduce the computational burden of the standard ITS method when closed-form expressions for the quantile function are unavailable, such as for custom distributions that incorporate assumptions about treatment effect waning as commonly encountered in practice,¹¹ and for flexible survival models, including relative survival models, spline models, mixture cure models, and response-based landmark models.^9,10 The implementation of our methods is identical across all possible survival models and can therefore be easily standardized from standard probabilistic decision analyses.

Generally, the precision of the EVSI estimator is influenced by the number of simulated data sets and the effective sample size of the simulated data. The discrete sampling method and, to a lesser degree, the interpolated ITS method, additionally introduce an approximation error that depends on the cycle length. It is generally recommended that discrete-time health economic models use a short cycle length to reduce the discrete-time approximation error, which could be as short as 1 week for slowly progressing chronic diseases.³¹ Our synthetic case study suggests that the approximation error introduced by our general-purpose methods is very small even when using a longer cycle length of 1 mo in combination with short follow-up times and a low effective sample size of the simulated data given a rapidly progressing disease.

We structured the synthetic case study around a PSM, a type of model that is frequently used to inform reimbursement decisions for new oncology drugs.²⁷ The key assumption behind a PSM is that survival endpoints, such as OS and PFS, are independent. This also implies that dependency between OS and PFS is not reflected in the EVSI data-simulation procedure when using a PSM. Joint modeling of OS and PFS could be implemented in a state transition model (STM), which uses transition probabilities to describe movements between health states over time. STMs require individual patient data to estimate all relevant transition probabilities, unlike PSMs, which can use digitized Kaplan–Meier data from published trials. OS and PFS data can be simulated jointly from a STM by first sampling a PFS time and then deciding whether the sampled PFS time is a progression or death event using a binomial experiment with probability derived from the hazards of transitioning from PFS to PPS and OS.³² If the sampled PFS time is a progression event, residual time until death can be simulated using the survival distribution for PPS to OS. Since transition probabilities are typically derived from survival curves fitted to time-to-event data, our data-simulation methods could also be useful in a STM framework.

If individual patient data are available with a similar study design and at least the same length of follow-up as the proposed study, study data sets could alternatively be simulated using a 2-level resampling method based on bootstrapping.³³ In this approach, the observed data set is first resampled $K$ times with replacement, and then $n$ values are sampled with replacement from the $k = 1, . . ., K$ resampled data sets, thereby generating $K$ new data sets with $n$ observations each. In most situations, particularly for novel treatments, this type of data will, however, not be available.

The key notion behind EVSI is that the prior distribution of the model parameters is updated with simulated study data to estimate the joint posterior distribution given both prior information and the simulated study data. In EVSI analyses, the prior distribution is often informed by external evidence, such as digitized Kaplan–Meier data. This may, however, not match the way in which real-world analyses of study data are conducted, since these may not synthesize the collected study data and the external evidence. Analysts should therefore ensure that the way in which the study data is analyzed once it has been collected is aligned with the assumptions underpinning the EVSI analysis.

Despite its routine application in cost-effectiveness analyses, there is currently a lack of guidance on how to model treatment effect waning.¹² In the synthetic case study, we modeled treatment effect waning by specifying probability distributions for the start and duration of the waning period, while preserving uncertainty about independent survival endpoints using a novel additive hazard approach. This had a large impact on the EVSI estimates, which highlights the importance of appropriately incorporating uncertainty about treatment effect waning in the EVSI calculations. There may be other possible approaches to model treatment effect waning, and these can easily be captured by our data-simulation methods as well.

Conclusion

The increasing prevalence of immature survival data in decision making, particularly for new cancer treatments,⁵ has been accompanied by the introduction of increasingly complex approaches for extrapolation,^9,10 which complicates the EVSI data-simulation step. Our general-purpose data-simulation methods greatly reduce the computational burden of the EVSI data-simulation step when custom distributions that incorporate treatment effect waning or flexible survival models are used for which closed-form expressions for the quantile function are unavailable. Our methods are straightforward to implement and can easily be automated from standard probabilistic decision analyses, such as those used in technology assessments of new pharmaceuticals.^11,34–36 This means that our general-purpose methods can be used to simulate survival data—with a similar accuracy and computational cost—as using the correct closed-form quantile function for any survival model. Efficient EVSI calculations for survival data can help decision makers determine whether current evidence is sufficient or whether there is a need for collecting additional survival data before making an adoption decision.

Supplemental Material

sj-docx-1-mdm-10.1177_0272989X231162069 – Supplemental material for General-Purpose Methods for Simulating Survival Data for Expected Value of Sample Information Calculations

Supplemental material, sj-docx-1-mdm-10.1177_0272989X231162069 for General-Purpose Methods for Simulating Survival Data for Expected Value of Sample Information Calculations by Mathyn Vervaart, Eline Aas, Karl P. Claxton, Mark Strong, Nicky J. Welton, Torbjørn Wisløff and Anna Heath in Medical Decision Making

Footnotes

The authors declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article. The authors disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: MV, EA, KC, MS, NJW, and TW were funded by a grant from the Norwegian Research Council through NordForsk (grant 298854). AH was supported by a Canada Research Chair in Statistical Trial Design and the Natural Sciences and Engineering Research Council of Canada (grant RGPIN-2021-03366) The funding agreement ensured the authors’ independence in designing the study, interpreting the data, writing, and publishing the report.

ORCID iDs

Mathyn Vervaart

Mark Strong

Torbjørn Wisløff

Anna Heath

Data Availability

The analysis code used in this study is available from Github at .

Supplemental Material

Supplementary material for this article is available online at .

References

Raiffa

Schlaifer

. Applied Statistical Decision Theory. Boston: Division of Research, Graduate School of Business Administration, Harvard University; 1961.

Vervaart

Strong

Claxton

Welton

Wisløff

Aas

. An efficient method for computing expected value of sample information for survival data from an ongoing trial. Med Decis Making. 2022;42:612–25.

Heath

Strong

Glynn

Kunst

Welton

Goldhaber-Fiebert

. Simulating study data to support expected value of sample information calculations: a tutorial. Med Decis Making. 2022;42:143–55.

Latimer

. Survival analysis for economic evaluations alongside clinical trials—extrapolation with patient-level data: inconsistencies, limitations, and a practical guide. Med Decis Making. 2013;33:743–54.

Tai

T-A

Latimer

Benedict

Kiss

Nikolaou

. Prevalence of immature survival data for anti-cancer drugs presented to the national institute for health and care excellence and impact on decision making. Value Health. 2021;24(4):505–12.

European Medicines Agency. Adaptive pathways. European Medicines Agency, 2018. Available from: https://www.ema.europa.eu/en/human-regulatory/research-development/adaptive-pathways. [Accessed 17 July , 2022].

European Medicines Agency. Conditional marketing authorisation. European Medicines Agency, 2018. Available from: https://www.ema.europa.eu/en/human-regulatory/marketing-authorisation/conditional-marketing-authorisation. [Accessed 17 July, 2022].

US Food and Drug Administration (FDA). Fast track, breakthrough therapy, accelerated approval, priority review. FDA, 2019. Available from: https://www.fda.gov/patients/learn-about-drug-and-device-approvals/fast-track-breakthrough-therapy-accelerated-approval-priority-review. [Accessed 17 July, 2022].

Rutherford

Lambert

Sweeting

, et al. NICE DSU Technical Support Document 21: flexible methods for survival analysis. 2020. p 97. Sheffield: Decision Support Unit, ScHARR, University of Sheffield. https://www.sheffield.ac.uk/nice-dsu/tsds/flexible-methods-survival-analysis

10.

Palmer

Borget

Friede

, et al. A guide to selecting flexible survival models to inform economic evaluations of cancer immunotherapies. Value Health. 2023;26(2):185–92.

11.

National Institute for Health and Care Excellence. Guide to the methods of technology appraisal 2013. 2013. Available from: https://www.nice.org.uk/process/pmg9

12.

Witlox

WJA

van Asselt

ADI

Wolff

, et al. Durvalumab for the treatment of locally advanced, unresectable, stage III non-small cell lung cancer: an evidence review group perspective of a NICE single technology appraisal. Pharmacoeconomics. 2020;38:317–24.

13.

Lambe

Duarte

Mahon

, et al. Cladribine tablets for the first-line treatment of relapsing-remitting multiple sclerosis: an evidence review group perspective of a NICE single technology Appraisal. Pharmacoeconomics. 2019;37:345–57.

14.

Strong

Oakley

Brennan

Breeze

. Estimating the expected value of sample information using the probabilistic sensitivity analysis sample: a fast, nonparametric regression-based method. Med Decis Making. 2015;35:570–83.

15.

Briggs

Claxton

Sculpher

. Decision Modelling for Health Economic Evaluation. 1st ed. Oxford: Oxford University Press; 2006.

16.

Ades

Claxton

. Expected value of sample information calculations in medical decision modeling. Med Decis Making. 2004;24:207–27.

17.

Collett

. Modelling Survival Data in Medical Research. 3rd ed. Boca Raton (FL): Chapman and Hall/CRC; 2015.

18.

R Core Team. R: A Language and Environment for Statistical Computing. R Foundation for Statistical Computing, Vienna, Austria, 2021. Available from: https://www.R-project.org/

19.

Jackson

Stevens

Ren

, et al. Extrapolating survival from randomized trials using external data: a review of methods. Med Decis Making. 2017;37:377–90.

20.

Kunst

Wilson

ECF

Glynn

, et al. computing the expected value of sample information efficiently: practical guidance and recommendations for four model-based methods. Value Health. 2020;23:734–42.

21.

Latimer

Adler

. Extrapolation beyond the end of trials to estimate long term survival and cost effectiveness. BMJ Med. 2022;1:1.

22.

Jackson

Thompson

Sharples

. Accounting for uncertainty in health economic decision models by using model averaging. J R Stat Soc Ser A Stat Soc. 2009;172:383–404.

23.

Jackson

Sharples

Thompson

. Structural and parameter uncertainty in Bayesian cost-effectiveness models. J R Stat Soc Ser C Appl Stat. 2010;59:233–53.

24.

Fragoso

Neto

. Bayesian model averaging: a systematic review and conceptual classiﬁcation. Int Stat Rev. 2018;86:1–28.

25.

Akaike

. Information theory and an extension of the maximum likelihood principle. In: Petrov

Csaki

, eds. Proceedings of the 2nd International Symposium on Information Theory. Budapest: Akademiai Kiado; 1973. p 267–81.

26.

Woods

Sideris

Palmer

Latimer

Soares

. Partitioned survival and state transition models for healthcare decision making in oncology: where are we now? Value Health. 2020;23:1613–21.

27.

Bojke

Soares

Claxton

, et al. Reference case methods for expert elicitation in health care decision making. Med Decis Making. 2022;42:182–93.

28.

Heath

Manolopoulou

Baio

. Efficient Monte Carlo estimation of the expected value of sample information using moment matching. Med Decis Making. 2018;38:163–73.

29.

Wood

. mgcv: Mixed GAM computation vehicle with automatic smoothness estimation, 2020. Available from: https://CRAN.R-project.org/package=mgcv

30.

Strong

Oakley

Brennan

. Estimating multiparameter partial expected value of perfect information from a probabilistic sensitivity analysis sample: a nonparametric regression approach. Med Decis Making. 2014;34:311–26.

31.

O’Mahony

Newall

van Rosmalen

. Dealing with time in health economic evaluation: methodological issues and recommendations for practice. Pharmacoeconomics. 2015;33:1255–68.

32.

Meller

Beyersmann

Ruﬁbach

. Joint modeling of progression-free and overall survival and computation of correlation measures. Stat Med. 2019;38:4270–89.

33.

Sadatsafavi

Marra

Bryan

. Two-level resampling as a novel method for the calculation of the expected value of sample information in economic trials. Health Econ. 2013;22:877–82.

34.

Canada’s Drug and Health Technology Agency (CADTH). Guidelines for the Economic evaluation of health technologies: Canada (4th edition), 2017. Available from: https://www.cadth.ca/sites/default/files/pdf/guidelines_for_the_economic_evaluation_of_health_technologies_canada_4th_ed.pdf

35.

Commonwealth of Australia. Guidelines for preparing a submission to the Pharmaceutical Beneﬁts Advisory Committee, version 5.0, 2016. Available from:https://pbac.pbs.gov.au/content/information/files/pbac-guidelines-version-5.pdf

36.

European Network for Health Technology Assessment. Practical considerations when critically as sessing economic evaluations, 2020. Available from: https://www.eunethta.eu/wp-content/uploads/2020/03/EUnetHTA-JA3WP6B2-5-Guidance-Critical-Assessment-EE_v1-0.pdf

Supplementary Material

Please find the following supplemental material available below.

For Open Access articles published under a Creative Commons License, all supplemental material carries the same license as the article it is associated with.

For non-Open Access articles published, all supplemental material carries a non-exclusive license, and permission requests for re-use of supplemental material or any part of supplemental material shall be sent directly to the copyright owner as specified in the copyright notice associated with the article.

0.00 MB

0.39 MB