The adaptr
package simulates adaptive clinical trials
using adaptive stopping, adaptive arm dropping and/or response-adaptive
randomisation.
The package has been developed as part of the INCEPT (Intensive Care Platform Trial) project, funded primarily by a grant from Sygeforsikringen “danmark”.
First, load the package:
library(adaptr)
Then, setup a trial with the desired specifications.
adaptr
offers the general purpose function
setup_trial()
, but here we use the built-in
setup_trial_binom()
for a trial with a binary, binomially
distributed, undesirable outcome such as mortality (adaptr
also includes setup_trial_norm()
for continuous, normally
distributed outcomes). The example trial specification has three
characteristics:
min_probs
).equivalence_prob
) of treatment differences
being < 5 %-points (equivalence_diff
), the trial is
stopped.soften_power
) by a
constant factor.<- setup_trial_binom(
binom_trial arms = c("Arm A", "Arm B", "Arm C"),
true_ys = c(0.25, 0.20, 0.30),
min_probs = rep(0.15, 3),
data_looks = seq(from = 300, to = 2000, by = 100),
equivalence_prob = 0.9,
equivalence_diff = 0.05,
soften_power = 0.5
)
See ?setup_trial()
for more details on these arguments
or vignette("Basic-examples", "adaptr")
for basic example
trial specifications and a thorough review of the general settings, and
vignette("Advanced-example", "adaptr")
for an
advanced example including details on how to specify
user-written functions for generating outcomes and posterior draws.
We can print an overview of the trial specification by simply running:
binom_trial#> Trial specification: generic binomially distributed outcome trial
#> * Undesirable outcome
#> * No common control arm
#> * Best arm: Arm B
#>
#> Arms, true outcomes, starting allocation probabilities
#> and allocation probability limits:
#> arms true_ys start_probs fixed_probs min_probs max_probs
#> Arm A 0.25 0.333 NA 0.15 NA
#> Arm B 0.20 0.333 NA 0.15 NA
#> Arm C 0.30 0.333 NA 0.15 NA
#>
#> Maximum sample size: 2000
#> Maximum number of data looks: 18
#> Planned data looks after: 300, 400, 500, 600, 700, 800, 900, 1000, 1100, 1200, 1300, 1400, 1500, 1600, 1700, 1800, 1900, 2000 patients
#>
#> Superiority threshold: 0.99
#> Inferiority threshold: 0.01
#> Equivalence threshold: 0.9 (no common control)
#> Absolute equivalence difference: 0.05
#> No futility threshold (not relevant - no common control)
#> Soften power for all analyses: 0.5
By default, probabilities are shown with 3 decimals. This can be
changed by explicitly print()
ing the specification with the
prob_digits
arguments, for example:
print(binom_trial, prob_digits = 2)
#> Trial specification: generic binomially distributed outcome trial
#> * Undesirable outcome
#> * No common control arm
#> * Best arm: Arm B
#>
#> Arms, true outcomes, starting allocation probabilities
#> and allocation probability limits:
#> arms true_ys start_probs fixed_probs min_probs max_probs
#> Arm A 0.25 0.33 NA 0.15 NA
#> Arm B 0.20 0.33 NA 0.15 NA
#> Arm C 0.30 0.33 NA 0.15 NA
#>
#> Maximum sample size: 2000
#> Maximum number of data looks: 18
#> Planned data looks after: 300, 400, 500, 600, 700, 800, 900, 1000, 1100, 1200, 1300, 1400, 1500, 1600, 1700, 1800, 1900, 2000 patients
#>
#> Superiority threshold: 0.99
#> Inferiority threshold: 0.01
#> Equivalence threshold: 0.9 (no common control)
#> Absolute equivalence difference: 0.05
#> No futility threshold (not relevant - no common control)
#> Soften power for all analyses: 0.5
Remember to define the seed
to ensure reproducible
results. The final trial results are wrapped here, but when used
interactively, you can widen the console for non-wrapped output:
<- run_trial(binom_trial, seed = 12345)
trial_res
trial_res#> Single simulation result: generic binomially distributed outcome trial
#> * Undesirable outcome
#> * No common control arm
#>
#> Final status: inconclusive, stopped at maximum sample size
#> Final/maximum allowed sample sizes: 2000/2000 (100.0%)
#>
#> Final trial results:
#> arms true_ys sum_ys ns raw_ests post_ests post_errs lo_cri hi_cri
#> Arm A 0.25 180 742 0.243 0.243 0.0159 0.213 0.275
#> Arm B 0.20 178 841 0.212 0.212 0.0141 0.185 0.241
#> Arm C 0.30 113 417 0.271 0.271 0.0215 0.230 0.316
#> final_status status_look status_probs final_alloc
#> active NA NA 0.194
#> active NA NA 0.656
#> inferior 2000 0.007 0.150
#>
#> Simulation details:
#> * Random seed: 12345
#> * Credible interval width: 95%
#> * Number of posterior draws: 5000
#> * Posterior estimation method: medians with MAD-SDs
Again, we can choose the number of decimals with
print()
:
print(trial_res, prob_digits = 2)
#> Single simulation result: generic binomially distributed outcome trial
#> * Undesirable outcome
#> * No common control arm
#>
#> Final status: inconclusive, stopped at maximum sample size
#> Final/maximum allowed sample sizes: 2000/2000 (100.0%)
#>
#> Final trial results:
#> arms true_ys sum_ys ns raw_ests post_ests post_errs lo_cri hi_cri
#> Arm A 0.25 180 742 0.24 0.24 0.016 0.21 0.27
#> Arm B 0.20 178 841 0.21 0.21 0.014 0.19 0.24
#> Arm C 0.30 113 417 0.27 0.27 0.022 0.23 0.32
#> final_status status_look status_probs final_alloc
#> active NA NA 0.19
#> active NA NA 0.66
#> inferior 2000 0.007 0.15
#>
#> Simulation details:
#> * Random seed: 12345
#> * Credible interval width: 95%
#> * Number of posterior draws: 5000
#> * Posterior estimation method: medians with MAD-SDs
Generally, we want to run many simulations using the same trial
specification to assess and compare performance metrics of different
trial designs. This is the job of run_trials()
(note the
final s); again, note the use of a
base_seed
for reproducible results. Here we run 25
simulations (for tolerable run-time when producing this example), but in
practice you will probably want to run more simulations.
<- run_trials(binom_trial, n_rep = 25, base_seed = 67890)
trial_res_mult
trial_res_mult#> Multiple simulation results: generic binomially distributed outcome trial
#> * Undesirable outcome
#> * Number of simulations: 25
#> * Number of simulations summarised: 25 (all trials)
#> * No common control arm
#> * Selection strategy: no selection if no superior arm
#> * Treatment effect compared to: no comparison
#>
#> Performance metrics (using posterior estimates):
#> * Sample sizes: mean 1324.0 (SD: 547.2) | median 1200.0 (IQR: 1000.0 to 1900.0)
#> * Total summarised outcomes: mean 296.1 (SD: 111.5) | median 282.0 (IQR: 228.0 to 415.0)
#> * Total summarised outcome rates: mean 0.228 (SD: 0.013) | median 0.230 (IQR: 0.220 to 0.235)
#> * Conclusive: 80.0%
#> * Superiority: 64.0%
#> * Equivalence: 16.0%
#> * Futility: 0.0% [not assessed]
#> * Inconclusive at max sample size: 20.0%
#> * Selection probabilities: Arm A: 0.0% | Arm B: 64.0% | Arm C: 0.0% | None: 36.0%
#> * RMSE: 0.01642
#> * RMSE treatment effect: not estimated
#> * Ideal design percentage: 100.0%
#>
#> Simulation details:
#> * Simulation time: 0.871 secs
#> * Base random seed: 67890
#> * Credible interval width: 95%
#> * Number of posterior draws: 5000
#> * Estimation method: posterior medians with MAD-SDs
run_trials()
can run simulations on several CPU cores
concurrently: set the cores
argument to some number greater
than 1
(which is the default value, resulting in
serial/non-parallel processing). As an aside, you can see the number of
available CPU cores by calling parallel::detectCores()
.
Just printing the resulting object (as the last line above) is not
necessarily what you want. adaptr
has two convenience
functions to get the results of multiple simulations of the same trial
specification, extract_results()
and
summary()
; see their documentation for details, e.g. on arm
selection in trials not ending with superiority. The
extract_results()
function yields a tidy
data.frame
with one simulation per row:
<- extract_results(trial_res_mult)
extr_res
nrow(extr_res)
#> [1] 25
head(extr_res)
#> sim final_n sum_ys ratio_ys final_status superior_arm selected_arm
#> 1 1 2000 415 0.2075000 max <NA> <NA>
#> 2 2 600 139 0.2316667 superiority Arm B Arm B
#> 3 3 1000 237 0.2370000 superiority Arm B Arm B
#> 4 4 900 209 0.2322222 equivalence <NA> <NA>
#> 5 5 2000 441 0.2205000 superiority Arm B Arm B
#> 6 6 1900 431 0.2268421 superiority Arm B Arm B
#> sq_err sq_err_te
#> 1 NA NA
#> 2 7.853843e-04 NA
#> 3 4.190319e-05 NA
#> 4 NA NA
#> 5 3.422824e-06 NA
#> 6 4.852161e-05 NA
summary()
(as expected and known from, e.g., regression
models in R
) summarises the results in a human-friendly
manner. Individual key results can be extracted from the resulting
object or printed:
<- summary(trial_res_mult)
res_sum
print(res_sum, digits = 1)
#> Multiple simulation results: generic binomially distributed outcome trial
#> * Undesirable outcome
#> * Number of simulations: 25
#> * Number of simulations summarised: 25 (all trials)
#> * No common control arm
#> * Selection strategy: no selection if no superior arm
#> * Treatment effect compared to: no comparison
#>
#> Performance metrics (using posterior estimates):
#> * Sample sizes: mean 1324.0 (SD: 547.2) | median 1200.0 (IQR: 1000.0 to 1900.0)
#> * Total summarised outcomes: mean 296.1 (SD: 111.5) | median 282.0 (IQR: 228.0 to 415.0)
#> * Total summarised outcome rates: mean 0.228 (SD: 0.013) | median 0.230 (IQR: 0.220 to 0.235)
#> * Conclusive: 80.0%
#> * Superiority: 64.0%
#> * Equivalence: 16.0%
#> * Futility: 0.0% [not assessed]
#> * Inconclusive at max sample size: 20.0%
#> * Selection probabilities: Arm A: 0.0% | Arm B: 64.0% | Arm C: 0.0% | None: 36.0%
#> * RMSE: 0.01642
#> * RMSE treatment effect: not estimated
#> * Ideal design percentage: 100.0%
#>
#> Simulation details:
#> * Simulation time: 0.871 secs
#> * Base random seed: 67890
#> * Credible interval width: 95%
#> * Number of posterior draws: 5000
#> * Estimation method: posterior medians with MAD-SDs
We need non-sparse results for plot_history()
(but
not for plot_status()
), so first we re-run
run_trials()
:
<- run_trials(binom_trial, n_rep = 25, base_seed = 67890,
trial_res_mult sparse = FALSE)
Note that the adaptr
plotting functions require the
ggplot2
package. First, we plot the overall trial statuses
according to the total number of patients randomised (this does
not require sparse = FALSE
):
plot_status(trial_res_mult, x_value = "total n")
We can also plot the statuses for specific arms. In this example, Arm A is increasingly dropped for inferiority as simulations progress:
plot_status(trial_res_mult, x_value = "total n", arm = "Arm A")
Next, we plot the history of allocation probabilities at each
adaptive analysis look. Intervals cover the inter-quartile range by
default (interval_width = 0.5
):
plot_history(trial_res_mult)
Plotting other summary metrics is possible; see
plot_history()
.
If using the package, please consider citing:
citation(package = "adaptr")
#>
#> To cite adaptr in publications use:
#>
#> Granholm A, Jensen AKG, Lange T, Kaas-Hansen BS (2022). adaptr: an R
#> package for simulating and comparing adaptive clinical trials.
#> Journal of Open Source Software, 7(72), 4284. URL
#> https://doi.org/10.21105/joss.04284.
#>
#> A BibTeX entry for LaTeX users is
#>
#> @Article{,
#> title = {{adaptr}: an R package for simulating and comparing adaptive clinical trials},
#> author = {Anders Granholm and Aksel Karl Georg Jensen and Theis Lange and Benjamin Skov Kaas-Hansen},
#> journal = {Journal of Open Source Software},
#> year = {2022},
#> volume = {7},
#> number = {72},
#> pages = {4284},
#> url = {https://doi.org/10.21105/joss.04284},
#> doi = {10.21105/joss.04284},
#> }