Skip to contents

Simulate dose-finding trials

sim_jointBLRM.Rd

Simulates dose-finding trials with up to six parallel monotherapy or two-drug combination therapy trials modeled together in a joint BLRM. The function assumes that two different compounds are involved, compounds 1 and 2 (and their combination). Up to two monotherapy trials for each compound can be actively simulated, and additionally two combination therapy trials. Note that the term "trial" can in this context also refer to trial arms (e.g. monotherapy and combination therapy), and that the underlying model will not differentiate between these notions. The order in which cohorts for the simulated trials are enrolled can be specified freely. Further, besides the actively simulated trials, the function supports the inclusion of historical data from arbitrarily many previous (already concluded) dose-finding trials for the considered compounds.

For an introduction to the use of this function, see the introduction vignette: vignette("intro_jointBLRM", package = "decider").

Note that the function sim_jointBLRM() uses the foreach framework to allow parallelized simulations. For this, a parallel backend needs to be registered by the user. A basic example of this is available in the following vignette: vignette("parallelization_jointBLRM", package = "decider").

The function expects a specification of dose levels and their assumed DLT rates for each involved trial, based on which the course of the dose-finding trial is simulated. Additional customization of the simulations is possible (separately for each simulated trial), e.g. regarding the escalation rule (EWOC and loss-based methods), definition of MTD, the maximum increment of escalations, the numbers of patients per cohort, and further specifications.

The six possible trials are encoded with the following suffixes: mono1.a, mono1.b, mono2.a, mono2.b, combi.a, and combi.b. The ending (.a or .b) differentiates between two different trials of the same kind, i.e., monotherapy for compound 1 (mono1), monotherapy for compound 2 (mono2), and combination therapy (combi). For some function arguments, one version of the argument is available for each trial (differentiated via their suffixes), so that different simulated trials may use different specifications in terms of e.g. dose levels, assumed DLT rate, and starting dose.

To differentiate between the up to six simulated trials and the actual simulation run that processes one realization of the activated trials, we will in the following use the term "study" for a single simulation run that determines a realization of each activated trial using the given specifications. For instance, if two trials are active, say mono1.a and combi.a, the term "trial" refers to the simulated trials with these names, while the notion "study" means a simulation run of both trials in the specified order and under the given input assumptions. Simulating e.g. 1000 studies means in this context that the trials mono1.a and combi.a are simulated in the specified order using the specified options for 1000 times.

Most of the parameters have sensible defaults and will often not need to be specified explicitly. The recommended workflow to set up the function is as follows:

  • Specify reference doses.

  • Optional: Specify available historical data and, if needed, specific priors for the model (the default ones are weakly informative).

  • Specify dose levels, starting doses and assumed DLT rates for each trial that needs to be simulated, and activate the corresponding trials via the active.[...] arguments.

  • Set the number of trials to be simulated.

  • Optional: Specify additional options for the simulations.

  • Optional: Register a parallel backend for the foreach framework, cf. foreach-package. Alternatively or additionally, when only a small number of cores is available, consider to supply a path for a working directory for saving and re-loading MCMC results. When a working directory is supplied, be aware that the function may delete or modify existing .RData files containing "_tmp" and the supplied file.name in their name from the working directory to avoid re-loading MCMC runs from previous simulations. Refer to the documentation below for more detail.

  • Run the simulation. Note that simulations of e.g. 1000 trials may take multiple hours depending on the specifications.

Usage

sim_jointBLRM(
   active.mono1.a = FALSE,
   active.mono1.b = FALSE,
   active.mono2.a = FALSE,
   active.mono2.b = FALSE,
   active.combi.a = FALSE,
   active.combi.b = FALSE,
   doses.mono1.a,
   doses.mono2.a,
   doses.mono1.b,
   doses.mono2.b,
   doses.combi.a,
   doses.combi.b,
   dose.ref1,
   dose.ref2,
   tox.mono1.a,
   tox.mono2.a,
   tox.combi.a,
   tox.mono1.b,
   tox.mono2.b,
   tox.combi.b,
   start.dose.mono1.a,
   start.dose.mono2.a,
   start.dose.mono1.b,
   start.dose.mono2.b,
   start.dose.combi.a1,
   start.dose.combi.a2,
   start.dose.combi.b1,
   start.dose.combi.b2,
   cohort.queue = rep( c(1,2,3,4,5,6), times= 100),
   historical.data = NULL,
   esc.rule = "ewoc",
   esc.comp.max=1,
   dosing.intervals = c(0.16, 0.33, 0.6),
   ewoc.threshold = 0.25,
   loss.weights = c(1, 0, 1, 2),
   dynamic.weights =  rbind(
       c(0.32, 0, 0.32, 0.36),
       c(0.29, 0, 0.31, 0.40),
       c(0.27, 0, 0.33, 0.40),
       c(0.20, 0, 0.30, 0.50)),
   prior.mu = list(
       mu_a1 = c(logit(0.33), 2),
       mu_b1 = c(0, 1),
       mu_a2 = c(logit(0.33), 2),
       mu_b2 = c(0, 1),
       mu_eta = c(0, 1.121)),
   prior.tau = list(
       tau_a1 = c(log(0.25), log(4)/1.96),
       tau_b1 = c(log(0.125), log(4)/1.96),
       tau_a2 = c(log(0.25), log(4)/1.96),
       tau_b2 = c(log(0.125), log(4)/1.96),
       tau_eta = c(log(0.125),log(4)/1.96)),
   saturating = FALSE,
   esc.step = NULL,
   esc.step.mono1.a = esc.step,
   esc.step.mono2.a = esc.step,
   esc.step.mono1.b = esc.step,
   esc.step.mono2.b = esc.step,
   esc.step.combi.a1 = esc.step,
   esc.step.combi.b1 = esc.step,
   esc.step.combi.a2 = esc.step,
   esc.step.combi.b2 = esc.step,
   esc.constrain = FALSE,
   esc.constrain.mono1.a=esc.constrain,
   esc.constrain.mono2.a=esc.constrain,
   esc.constrain.mono1.b=esc.constrain,
   esc.constrain.mono2.b=esc.constrain,
   esc.constrain.combi.a1=esc.constrain,
   esc.constrain.combi.b1=esc.constrain,
   esc.constrain.combi.a2=esc.constrain,
   esc.constrain.combi.b2=esc.constrain,
   cohort.size = c(3),
   cohort.size.mono1.a = cohort.size,
   cohort.size.mono1.b = cohort.size,
   cohort.size.mono2.a = cohort.size,
   cohort.size.mono2.b = cohort.size,
   cohort.size.combi.a = cohort.size,
   cohort.size.combi.b = cohort.size,
   cohort.prob = NULL,
   cohort.prob.mono1.a = cohort.prob,
   cohort.prob.mono1.b = cohort.prob,
   cohort.prob.mono2.a = cohort.prob,
   cohort.prob.mono2.b = cohort.prob,
   cohort.prob.combi.a = cohort.prob,
   cohort.prob.combi.b = cohort.prob,
   max.n = 42,
   max.n.mono1.a = max.n,
   max.n.mono1.b = max.n,
   max.n.mono2.a = max.n,
   max.n.mono2.b = max.n,
   max.n.combi.a = max.n,
   max.n.combi.b = max.n,
   mtd.decision = list(
       target.prob = 0.5,
       pat.at.mtd = 6,
       min.pat = 12,
       min.dlt = 1,
       rule = 2),
   mtd.decision.combi.a = mtd.decision,
   mtd.decision.combi.b = mtd.decision,
   mtd.decision.mono1.a = mtd.decision,
   mtd.decision.mono1.b = mtd.decision,
   mtd.decision.mono2.a = mtd.decision,
   mtd.decision.mono2.b = mtd.decision,
   mtd.enforce = FALSE,
   mtd.enforce.mono1.a = mtd.enforce,
   mtd.enforce.mono2.a = mtd.enforce,
   mtd.enforce.mono1.b = mtd.enforce,
   mtd.enforce.mono2.b = mtd.enforce,
   mtd.enforce.combi.a = mtd.enforce,
   mtd.enforce.combi.b = mtd.enforce,
   backfill.mono1.a = FALSE,
   backfill.mono1.b = FALSE,
   backfill.mono2.a = FALSE,
   backfill.mono2.b = FALSE,
   backfill.combi.a = FALSE,
   backfill.combi.b = FALSE,
   backfill.size = c(3),
   backfill.prob = NULL,
   backfill.size.mono1.a = backfill.size,
   backfill.size.mono1.b = backfill.size,
   backfill.size.mono2.a = backfill.size,
   backfill.size.mono2.b = backfill.size,
   backfill.size.combi.a = backfill.size,
   backfill.size.combi.b = backfill.size,
   backfill.prob.mono1.a = backfill.prob,
   backfill.prob.mono1.b = backfill.prob,
   backfill.prob.mono2.a = backfill.prob,
   backfill.prob.mono2.b = backfill.prob,
   backfill.prob.combi.a = backfill.prob,
   backfill.prob.combi.b = backfill.prob,
   backfill.start.mono1.a = NULL,
   backfill.start.mono1.b = NULL,
   backfill.start.mono2.a = NULL,
   backfill.start.mono2.b = NULL,
   backfill.start.combi.a1 = NULL,
   backfill.start.combi.a2 = NULL,
   backfill.start.combi.b1 = NULL,
   backfill.start.combi.b2 = NULL,
   n.studies = 1,
   seed = sample.int(.Machine$integer.max, 1),
   chains = 4,
   iter = 13500,
   warmup = 1000,
   adapt_delta = 0.8,
   max_treedepth = 15,
   refresh=0,
   file.name = NULL,
   path = NULL,
   monitor.path = NULL,
   working.path = NULL,
   clean.working.path = FALSE,
   output.sim.config = TRUE,
   output.sim.detail = FALSE
)

Arguments

active.mono1.a, active.mono1.b, active.mono2.a, active.mono2.b, active.combi.a, active.combi.b

Logicals, default to FALSE. Each parameter activates the simulation of the corresponding trial. If the parameter is FALSE, the corresponding trial is not simulated, all parameters with this suffix are ignored, and no historical information can be included for this trial.

doses.mono1.a, doses.mono1.b, doses.mono2.a, doses.mono2.b

Numericals, one dimensional vectors with positive, strictly ascending entries. Each vector indicates the dose levels for the corresponding monotherapy trial. Optional if active.[...] is FALSE for the corresponding suffix.

doses.combi.a, doses.combi.b

Numericals, two dimensional arrays with two rows, arbitrarily many columns, and positive entries. Does not need to be ascending. Each column indicates a combination dose level for the corresponding combination therapy trial, where the entry in the first row of a specific column gives the dose level of compound 1, and the entry in the second row the dose level of compound 2. Optional if active.combi.[...] is FALSE for the corresponding suffix.

dose.ref1

Positive numerical value. Reference dose of compound 1.

dose.ref2

Positive numerical value. Reference dose of compound 2.

tox.mono1.a, tox.mono1.b, tox.mono2.a, tox.mono2.b

Numericals, one dimensional vectors with entries in \((0,1)\). Optional if active.mono.[...] is FALSE for the corresponding suffix. Otherwise, must have the same length as doses.[...] for the corresponding suffix. The entry \(i\) specifies the DLT rate of dose \(i\) of the corresponding trial.

tox.combi.a, tox.combi.b

Numericals, one dimensional vectors with entries in \((0,1)\). Optional if active.combi.[...] is FALSE for the corresponding suffix. Otherwise, its length must equal the number of columns of doses.combi.[...]. The entry \(i\) specifies the DLT rate of the dose combination of row \(i\) of doses.combi.[...].

start.dose.mono1.a, start.dose.mono2.a, start.dose.mono1.b, start.dose.mono2.b

Positive numerical values that specify the starting dose for the simulated monotherapy trials. The value of start.dose.mono.[...] must be one of the entries of doses.mono.[...] for the corresponding suffix. The starting doses are optional if the corresponding trial is not activated.

start.dose.combi.a1, start.dose.combi.a2, start.dose.combi.b1, start.dose.combi.b2

Positive numerical values that specify the starting dose for the simulated combination therapy trials. The values with suffixes .a1 and .a2 define the starting dose for first and second compound of the combi.a trial, and must appear as combination in one of the columns of doses.combi.a. The same applies for the trial with suffix combi.b. The starting doses are optional if the corresponding trial is not activated.

cohort.queue

Optional numerical or character vector that specifies the order or pattern in which cohorts are enrolled in the simulated trials. By default, cohorts are enrolled sequentially in the order: mono1.a, mono2.a, combi.a, mono1.b, mono2.b, combi.b, which is repeated until all trials concluded. In general, entry \(i\) specifies to which simulated trial the \(i\)th simulated cohort is enrolled (provided that this trial has not concluded yet, in which case the entry would be ignored). The cohort.queue can therefore be used to express e.g. that a certain number of monotherapy cohorts are completed prior to the start of a combination therapy trial. The cohort.queue can either be a short pattern that is repeated automatically, or alternatively a long vector that contains already the order of all cohorts that may occur in the trial. The latter can e.g. be used to ensure that a trial has concluded prior to the start of some other trial by supplying enough cohorts to the former to reach its maximum sample size before the first cohort for the latter trial is enrolled.

The trials to which cohorts are assigned are encoded as follows:

  • "mono1.a" or 1
    Corresponds to a cohort for the mono1.a trial (monotherapy trial for compound 1).

  • "mono2.a" or 2
    Corresponds to a cohort for the mono2.a trial (monotherapy trial for compound 2).

  • "combi.a" or 3
    Corresponds to a cohort for the combi.a trial (combination therapy trial).

  • "mono1.b" or 4
    Corresponds to a cohort for the mono1.b trial (monotherapy trial for compound 1).

  • "mono2.b" or 5
    Corresponds to a cohort for the mono2.b trial (monotherapy trial for compound 2).

  • "combi.b" or 6
    Corresponds to a cohort for the combi.b trial (combination therapy trial).

Assuming that each cohort is assigned to the minimal possible cohort size (specified via the argument cohort.size), the cohort.queue should contain sufficiently many entries for each activated trial so that it can in principle reach the specified maximum patient number to be enrolled in the trial. For instance, if cohorts consist of 3 or more patients and the maximum sample size (given in max.n.[...]) of this trial is 30, than the trial should appear at least 10 times in the cohort queue. If this is not the case, the function will repeat the cohort queue until sufficiently many cohorts are contained.

historical.data

Optional parameter that must be NULL (the default) or a named list that specifies the historical data to be included in the trial. Historical data can be included from both purely historical trials (i.e., not actively simulated) and also directly to one of the simulated trials (in which case the simulation assumes that this data was previously recorded in the corresponding trial).

A value of NULL represents that no historical data is available. Otherwise, the historical data is specified via the named list entries dose1, dose2, n.pat, n.dlt, and trial (where the names are not case-sensitive). The list entries should be vectors whose length is equal to the number of historical cohorts to be included. Their format should be as follows:

  • historical.data$dose1
    Numerical vector with non-negative or NA entries. Entry \(i\) gives the dose level of compound 1 that was administered to the \(i\)th historical cohort. If compound 1 was not administered to this cohort, the entry should be 0 or NA.

  • historical.data$dose2
    Numerical vector with non-negative or NA entries. Entry \(i\) gives the dose level of compound 2 that was administered to the \(i\)th historical cohort. If compound 2 was not administered to this cohort, the entry should be 0 or NA.

  • historical.data$n.pat
    Numerical vector with non-negative entries. Entry \(i\) gives the number of patients of the \(i\)th historical cohort. If the patient number is 0, the historical cohort is removed automatically. NA entries are not permitted.

  • historical.data$n.dlt
    Numerical vector with non-negative entries. Entry \(i\) gives the number of patients that experienced DLT among the \(i\)th historical cohort. The number of DLTs is not allowed be larger than the number of patients in the cohort. NA entries are not permitted.

  • historical.data$trial
    Numerical or character vector with arbitrary entries (except for NA which is not permitted). Entry \(i\) gives the code for the trial from which the \(i\)th historical cohort arises. That is, historical patients are grouped in different historical or concurrent trials according to these values, where the entries are not case sensitive and do not differentiate between numbers and the converted strings (i.e., 42 and "42" represent the same trial). The following entries are reserved to represent that the cohort belongs to one of the actively simulated trials:

    • "mono1.a" or 1
      The cohort is assumed to belong to the actively simulated trial mono1.a.

    • "mono2.a" or 2
      The cohort is assumed to belong to the actively simulated trial mono2.a.

    • "combi.a" or 3
      The cohort is assumed to belong to the actively simulated trial combi.a.

    • "mono1.b" or 4
      The cohort is assumed to belong to the actively simulated trial mono1.b.

    • "mono2.b" or 5
      The cohort is assumed to belong to the actively simulated trial mono2.b.

    • "combi.b" or 6
      The cohort is assumed to belong to the actively simulated trial combi.b.

    All other strings or numbers are interpreted to belong to purely historical (already completed) trials. If data for the actively simulated trials are included, the corresponding trial must be activated (otherwise, the reserved trial codes are not permitted), and the historical cohorts must match the type of the trial. The latter means that historical data for e.g. the mono1.a trial must contain observations from monotherapy with compound 1.
esc.rule

Optional character string, must have one of the following values: "ewoc", "ewoc.opt", "ewoc.max", "loss", "dynamic", or "dynamic.loss", where the default is "ewoc". The value of esc.rule specifies which escalation rule is applied during the simulation, cf. Details. The character strings are interpreted in a case-insensitive fashion. The possible values are interpreted as follows:

  • "ewoc" or "ewoc.opt"
    Applies the EWOC criterion with optimal probability rule and three dosing intervals (underdosing, target dosing, and overdosing). That is, the next dose is the one with the largest probability that the true DLT rate lies in the target interval among the doses whose probability of having a DLT rate in the overdosing interval is lower than the feasibility bound specified in the argument ewoc.threshold. If multiple doses or dose combinations attain the maximum target probability and satisfy the EWOC principle, the escalation rule will additionally maximize the dose level of the compound that is specified in the argument esc.comp.max.

  • "ewoc.max"
    Applies the EWOC criterion with maximum dose rule and three dosing intervals (underdosing, target dosing, and overdosing). That is, the next dose is the largest one among the doses whose probability of having true DLT rate in the overdosing interval is lower than the feasibility bound specified in the argument ewoc.threshold. For dose combinations, the escalation rule will first maximize the compound specified in esc.comp.max, and afterwards the other compound.

  • "loss"
    Applies the static loss escalation rule using the weights specified in the argument loss.weights and four dosing intervals (under dosing, target dosing, excessive toxicity, and unacceptable toxicity). When multiple dose combinations result in the same expected loss (Bayes risk), the function will additionally maximize the compound specified in the argument esc.comp.max among the doses with minimum expected loss.

  • "dynamic" or "dynamic.loss"
    Applies the dynamic loss escalation rule using the weights specified in the argument dynamic.weights and four dosing intervals (under dosing, target dosing, excessive toxicity, and unacceptable toxicity). When multiple dose combinations result in the same expected loss (Bayes risk), the function will additionally maximize the compound specified in the argument esc.comp.max among the doses with minimum expected loss.

esc.comp.max

Optional integer, must be either 1 (the default) or 2. The value specifies which compound (1 or 2) is maximized first when an escalation rule leads to a draw among multiple dose combinations.

dosing.intervals

Optional numeric vector with ascending entries between 0 and 1. The argument dosing.intervals must be of length 1, 2, or 3, and defines the boundaries of the dosing intervals used for escalation decisions. The default is c(0.16, 0.33, 0.6) which can be used for any escalation rule. The interpretation and requirements vary slightly across different escalation rules (different values of esc.rule).

More precisely, the following cases are differentiated. If esc.rule is...

  • "ewoc","ewoc.max", or "ewoc.opt"
    The length can be either 1, 2, or 3. In the latter two cases, only the first two entries are used to specify the targeted dosing interval. That is, the target interval is assumed to range from dosing.intervals[1] to dosing.intervals[2], while a potential third entry is discarded. If dosing.intervals has length 1 and value x, the argument is internally transformed to c(0, x). That is, the target interval is assumed to range from 0 to x, while the underdosing interval is not considered.
    The overdosing interval is implicitly defined to range from dosing.intervals[2] to 1, while the underdosing interval is defined to range from 0 to dosing.intervals[1].

  • "loss","dynamic.loss", or "dynamic"
    Must have three entries which are used to specify the targeted dosing interval and the interval for excessively toxic dosing. That is, the target interval is assumed to range from dosing.intervals[1] to dosing.intervals[2], while the interval for excessively toxic dosing is assumed to range from dosing.intervals[2] to dosing.intervals[3]. The underdosing interval is therefore implicitly defined to range from 0 to dosing.intervals[1], while the unacceptably toxic dosing interval is implicitly defined to range from dosing.intervals[3] to 1.

ewoc.threshold

Optional numerical value between 0 and 1 (excluding the boundaries), defaults to 0.25. Defines the feasibility bound for the EWOC criterion, i.e., when esc.rule is one of the EWOC-based rules, the next dose must have a probability of having DLT rate in the overdosing interval lower than ewoc.threshold. The value is ignored when esc.rule specifies one of the loss-based escalation rules.

loss.weights

Optional numerical vector with four entries (which can be arbitrary numbers), the default is c(1,0,1,2). Specifies the interval weights/penalties that are used for static loss escalation. This is only needed when esc.rule is "loss", and the parameter is ignored otherwise. More precisely, loss.weights[1] is the weight of the underdosing interval, loss.weights[2] the weight of the target dosing interval, loss.weights[3] the weight of the excessively toxic dosing interval, and loss.weights[4] is the weight of the unacceptably toxic dosing interval. See the Details section below for explanations regarding the methodology.

dynamic.weights

Optional numerical matrix with four rows and four columns, and arbitrary numbers as entries. Specifies the interval weights/penalties that are used for dynamic loss escalation. This is only needed when esc.rule is "dynamic" or "dynamic.loss", and the parameter is ignored otherwise. Each row of dynamic.weights specifies one of the static loss weight vectors that are interpolated during dynamic loss escalation. See the Details section below for explanations regarding the methodology.

More precisely, dynamic.weights[1, ] is the static weight vector that is weighted by the posterior probability that the reference dose has true DLT rate in the underdose interval, dynamic.weights[2, ] is the static weight vector that is weighted by the posterior probability that the reference dose has true DLT rate in the target interval, dynamic.weights[3, ] is the static weight vector that is weighted by the posterior probability that the reference dose has true DLT rate in the excessively toxic dosing interval, and dynamic.weights[4, ] is the static weight vector that is weighted by the posterior probability that the reference dose has true DLT rate in the unacceptably toxic dosing interval. The default value for dynamic.weights is a matrix with the following four rows (in order):

  • c(0.32, 0, 0.32, 0.36),

  • c(0.29, 0, 0.31, 0.40),

  • c(0.27, 0, 0.33, 0.40),

  • c(0.20, 0, 0.30, 0.50).

prior.mu

Optional list that gives the prior distribution for the hyper means \(\mu\). The list must have named entries, which all need to be numeric vectors of length 2:

  • prior.mu$mu_a1
    Numeric with length two, defaults to c(logit(0.33), 2). Specifies mean and SD of the hypermean \(\mu_1\) of the parameter \(log(\alpha_1)\) in the BLRM. The second entry must therefore be positive. See the Details section below for more detail.

  • prior.mu$mu_b1
    Numeric with length two, defaults to c(0, 1). Specifies mean and SD of the hypermean \(\mu_2\) of the parameter \(log(\beta_1)\) in the BLRM. The second entry must therefore be positive. See the Details section below for more detail.

  • prior.mu$mu_a2
    Numeric with length two, defaults to c(logit(0.33), 2). Specifies mean and SD of the hypermean \(\mu_3\) of the parameter \(log(\alpha_2)\) in the BLRM. The second entry must therefore be positive. See the Details section below for more detail.

  • prior.mu$mu_b2
    Numeric with length two, defaults to c(0, 1). Specifies mean and SD of the hypermean \(\mu_4\) of the parameter \(log(\beta_2)\) in the BLRM. The second entry must therefore be positive. See the Details section below for more detail.

  • prior.mu$mu_eta
    Numeric with length two, defaults to c(0, 1.121). Specifies mean and SD of the hypermean \(\mu_5\) of the parameter \(\eta\) in the BLRM. The second entry must therefore be positive. See the Details section below for more detail.

prior.tau

Optional list that gives the prior distribution for the between-trial heterogeneities (hyper SD) \(\tau\). The list must have named entries, which all need to be numeric vectors of length 2:

  • prior.tau$tau_a1
    Numeric with length two, defaults to c(log(0.25), log(2)/1.96). Specifies mean and SD on log-scale of the heterogeneity \(\tau_1\) of the parameter \(log(\alpha_1)\) in the BLRM. The second entry must therefore be positive. See the Details section below for more detail.

  • prior.tau$tau_b1
    Numeric with length two, defaults to c(log(0.125), log(2)/1.96). Specifies mean and SD on log-scale of the heterogeneity \(\tau_2\) of the parameter \(log(\beta_1)\) in the BLRM. The second entry must therefore be positive. See the Details section below for more detail.

  • prior.tau$tau_a2
    Numeric with length two, defaults to c(log(0.25), log(2)/1.96). Specifies mean and SD on log-scale of the heterogeneity \(\tau_3\) of the parameter \(log(\alpha_2)\) in the BLRM. The second entry must therefore be positive. See the Details section below for more detail.

  • prior.tau$tau_b2
    Numeric with length two, defaults to c(log(0.125), log(2)/1.96). Specifies mean and SD on log-scale of the heterogeneity \(\tau_4\) of the parameter \(log(\beta_2)\) in the BLRM. The second entry must therefore be positive. See the Details section below for more detail.

  • prior.tau$tau_eta
    Numeric with length two, defaults to c(log(0.125), log(2)/1.96). Specifies mean and SD on log-scale of the heterogeneity \(\tau_5\) of the parameter \(\eta\) in the BLRM. The second entry must therefore be positive. See the Details section below for more detail.

saturating

Optional logical, defaults to FALSE. If TRUE, the BLRM will be using a saturating interaction term as described in OncoBayes2::blrm_formula_saturating(). Also refer to the Details section in the documentation of scenario_jointBLRM().

esc.step, esc.step.mono1.a, esc.step.mono2.a, esc.step.combi.a1, esc.step.combi.a2

Optional numerical values that specify the maximum factor of dose escalations that is demanded additionally to the selected escalation rule. The default is NULL, in which case the escalation step is determined automatically for the corresponding trial (if activated) as the maximum factor between the (sorted) available doses for the trial. As an example, an value of esc.step=2 indicates that the dose can at most be doubled, or, in other words, that the increment between subsequent dose levels is at most 100%.

The argument esc.step is a global variant, all escalation steps that were not specified by the user default to its value. In particular, if all trials shall use the same escalation factor, only the value of esc.step needs to be defined, while the trial-specific variants can be used to express that a specific trial shall deviate from the global parameter. Note that consistency checks are performed for these parameters, i.e., whether the specified escalation step allows to reach all available doses.

The value of esc.step.[...] gives the escalation step that is used for the corresponding trial. Note that there are two escalation steps for combination therapy trials: For instance, esc.step.combi.a1 gives the maximum escalation factor in the first compound during combination therapy escalations, while esc.step.combi.a2 defines the escalation step for the second compound.

esc.step.mono1.b, esc.step.mono2.b, esc.step.combi.b1, esc.step.combi.b2

Same as esc.step.[...].a (where [...] is mono1, mono2, or combi) but for the second potentially simulated trial (suffix .b) of the respective trial type.

esc.constrain, esc.constrain.mono1.a, esc.constrain.mono2.a

Optional logicals, default to FALSE. When TRUE, escalation decisions in the corresponding trial are constrained to the available dose levels, in the sense that only the next larger dose can be used for the next cohort. This does not apply to de-escalation, where arbitrary lower dose levels are allowed. Note that esc.constrain[...] overrules esc.step[...], so that the given escalation step is ignored when esc.constrain[...] is TRUE.

As before, esc.contrain is the global variant and affects all trials. If a specific trial shall deviate from the global setting, this needs to be activated specifically with the trial specific variant esc.constrain.[...] using the corresponding suffix.

esc.constrain.mono1.b, esc.constrain.mono2.b, esc.constrain.combi.b1, esc.constrain.combi.b2

Same as esc.constrain.[...].a (where [...] is mono1, mono2, or combi) but for the second potentially simulated trial (suffix .b) of the respective trial type.

esc.constrain.combi.a1, esc.constrain.combi.a2

Optional logicals, default to FALSE. Work in the same way as esc.constrain.mono[...], except that the compounds in a combination therapy trial are constrained separately from each other, where e.g. esc.constrain.combi.a1 affects compound 1, while esc.constrain.combi.a2 affects compound 2. For instance, it is possible to activate the restriction to the next larger dose level only for the first compound, and still use an escalation step for the second compound.

cohort.size, cohort.size.mono1.a, cohort.size.mono2.a, cohort.size.combi.a

Optional positive integer vectors that specify the available cohort sizes for simulated cohorts. The global variant, cohort.size, defaults to 3 (i.e., all cohorts will contain 3 patients), the remaining trial-specific parameters default to the value of cohort.size. In particular, if all trials shall use the same cohort sizes, only the value of cohort.size needs to be defined, while the trial-specific variants can be used to express that a specific trial shall deviate from the global parameter.

The cohort.size can either be a single value (all cohorts use the same size) or a vector of integers that give the available cohort sizes. In the latter case, the simulated cohort sizes are drawn at random from the possible cohort sizes using the probability for each cohort size that is specified in cohort.prob. As before, the trial-specific variant, cohort.size.[...], only affect the corresponding trial.

cohort.size.mono1.b, cohort.size.mono2.b, cohort.size.combi.b

Same as cohort.size.[...].a (where [...] is mono1, mono2, or combi) but for the second potentially simulated trial (suffix .b) of the respective trial type.

cohort.prob, cohort.prob.mono1.a, cohort.prob.mono2.a, cohort.prob.combi.a

Optional positive numeric vectors with values between 0 and 1 that specify the probability for each of the available cohort sizes in the corresponding argument cohort.size.[...]. The global variant, cohort.prob, defaults to NULL, the remaining trial-specific parameters default to the value of cohort.prob. If the parameter is NULL for a trial, the function will assign the same probability to each entry of cohort.size, i.e., the cohort size is chosen uniformly at random from the provided possible cohort sizes. If all trials shall use the same probabilities for selecting the cohort size, only the value of cohort.prob needs to be defined, while the trial-specific variants can be used to express that a specific trial shall deviate from the global parameter.

Each of the arguments cohort.prob.[...] must be NULL or a vector of the same length as the argument cohort.size.[...], where the lengths and contents may differ across different trials/suffixes.

cohort.prob.mono1.b, cohort.prob.mono2.b, cohort.prob.combi.b

Same as cohort.prob.[...].a (where [...] is mono1, mono2, or combi) but for the second potentially simulated trial (suffix .b) of the respective trial type.

max.n, max.n.mono1.a, max.n.mono2.a, max.n.combi.a

Optional positive integer values that specify the maximum number of patients to be enrolled in simulated trials. If the simulated number of patients in a trial equals or surpasses the specified number in max.n.[...] and has not stopped prematurely and not determined an MTD yet, the simulated trial is stopped without declaring an MTD due to reaching the maximum sample size.

The global variant of the parameter, max.n, defaults to 42 (i.e., MTD must be reached after 14 cohorts when the default values for cohort.size and cohort.prob are used). The trial-specific parameters, max.n.[...], default to the global value max.n and only need to be adjusted if the corresponding trial shall deviate from the global value.

max.n.mono1.b, max.n.mono2.b, max.n.combi.b

Same as max.n.[...].a (where [...] is mono1, mono2, or combi) but for the second potentially simulated trial (suffix .b) of the respective trial type.

mtd.decision, mtd.decision.mono1.a, mtd.decision.mono2.a, mtd.decision.combi.a

Optional named lists that specify the rules for MTD selection. The determination of the MTD may occur independently for each simulated trial as soon as the conditions that are specified in the corresponding mtd.decision[...] are met. Each list must have five named entries, namely, target.prob, pat.at.mtd, min.pat, min.dlt, and rule. The names of the list entries are not case sensitive. For the actual MTD selection, the function will always first demand that the current dose level of the current trial is recommended again, i.e. is also the recommended dose for the next cohort in this trial (the so-called stabilization criterion). If this is the case, the conditions specified in the corresponding mtd.decision[...] are checked. If these conditions are also satisfied, the function will declare the current dose of this trial to be the MTD and will not simulate any further cohorts for this trial.

The entries are single numerical values with the following specifications and interpretations.

  • mtd.decision[...]$rule
    Integer, can only be 1 or 2 (the latter is the default). Specifies which of the remaining entries of mtd.decision[...] are used as mandatory conditions for MTD selection. For both rules, two conditions are demanded as necessary conditions, while only one of the remaining two conditions must hold. The precise specifications are:

    • mtd.decision[...]$rule=1
      The conditions posed by the entries mtd.decision[...]$pat.at.mtd and mtd.decision[...]$min.dlt are mandatory, while only one of the conditions that are posed by mtd.decision[...]$min.pat and mtd.decision[...]$target.prob must hold.

    • mtd.decision[...]$rule=2
      The conditions from mtd.decision[...]$min.pat and mtd.decision[...]$min.dlt are mandatory, while only one of the conditions that are posed by mtd.decision[...]$pat.at.mtd and mtd.decision[...]$target.prob must hold.

  • mtd.decision[...]$min.pat
    Non-negative integer, defaults to 12. Specifies the minimum number of patients that need to have been enrolled in a trial before the MTD can be declared.

  • mtd.decision[...]$pat.at.mtd
    Non-negative integer, defaults to 6. Specifies the minimum number of patients that need to have been treated with the current dose (i.e., the current candidate for MTD determination) before the MTD can be declared.

  • mtd.decision[...]$min.dlt
    Non-negative integer, defaults to 1. Specifies the minimum number of DLTs that need to have been observed in a trial before the MTD can be declared.

  • mtd.decision[...]$target.prob
    Non-negative real number strictly smaller than 1, defaults to 0.5. Specifies a minimum probability of having true DLT rate in the target interval that needs to be reached by the dose that is to be declared the MTD. Note that this condition is always one of the two optional conditions, i.e., in case of the first rule, the condition does not need to hold when sufficiently many patients were enrolled in total, while in case of the second rule, the condition does not need to hold when sufficiently many patients were enrolled at the MTD.

As before, the global parameter mtd.decision can be used to select a rule that is used for all simulated trials. The trial-specific parameters, mtd.decision.[...], default to the global value and can be adjusted to use different, trial-specific decision rules for some or all of the trials.

mtd.decision.mono1.b, mtd.decision.mono2.b, mtd.decision.combi.b

Same as mtd.decision.[...].a (where [...] is mono1, mono2, or combi) but for the second potentially simulated trial (suffix .b) of the respective trial type.

mtd.enforce, mtd.enforce.mono1.a, mtd.enforce.mono2.a, mtd.enforce.combi.a

Optional logicals, default to FALSE. These parameters indicate whether the declaration of the MTD should be enforced once the maximum number of patients is reached and the trial was not stopped yet. That is, if the parameter mtd.enforce.[...] for a specific trial is TRUE, the function will declare the next recommended dose as MTD once the maximum number of patients (max.n.[...]) was reached (provided that the trial was not stopped early or found the MTD before reaching the maximum sample size). In particular, these trials are not counted to belong to the category of trials that reached the maximum sample size without declaring an MTD. Further, note that the value of mtd.decision.[...]$min can be set to a value greater than the maximum sample size, which would result in only declaring an MTD if the maximum sample size is reached, but not at earlier points (as the condition on the minimal patient number could not be fulfilled then).

As before, there is a global variant of this parameter, namely mtd.enforce, which can be used to set this option for all trials at once. To set the escalation rule only for a specific trial, use the trial-specific variants mtd.enforce.[...] with the corresponding suffix.

mtd.enforce.mono1.b, mtd.enforce.mono2.b, mtd.enforce.combi.b

Same as mtd.enforce.[...].a (where [...] is mono1, mono2, or combi) but for the second potentially simulated trial (suffix .b) of the respective trial type.

backfill.mono1.a, backfill.mono2.a, backfill.combi.a

Optional logicals, indicating whether back-fill cohorts are simulated for the trial.

backfill.mono1.b, backfill.mono2.b, backfill.combi.b

Same as backfill.[...].a (where [...] is mono1, mono2, or combi) but for the second potentially simulated trial (suffix .b) of the respective trial type.

backfill.size, backfill.size.mono1.a, backfill.size.mono2.a, backfill.size.combi.a

Optional numericals, provide the size of simulated back-fill cohorts. Interpreted in the same fashion as cohort.size.

backfill.prob, backfill.prob.mono1.a, backfill.prob.mono2.a, backfill.prob.combi.a

Optional numericals, provide the probabilities if multiple possible back-fill cohort sizes are given. Interpreted in the same fashion as cohort.prob.

backfill.size.mono1.b, backfill.size.mono2.b, backfill.size.combi.b

Same as backfill.size.[...].a (where [...] is mono1, mono2, or combi) but for the second potentially simulated trial (suffix .b) of the respective trial type.

backfill.prob.mono1.b, backfill.prob.mono2.b, backfill.prob.combi.b

Same as backfill.prob.[...].a (where [...] is mono1, mono2, or combi) but for the second potentially simulated trial (suffix .b) of the respective trial type.

backfill.start.mono1.a, backfill.start.mono2.a, backfill.start.combi.a1, backfill.start.combi.a2

Optional numericals. Specify the first dose on which back-fill cohorts are enrolled. If not provided, the lowest available dose will be assumed to be the starting point of back-fill cohorts in the respective trial.

backfill.start.mono1.b, backfill.start.mono2.b, backfill.start.combi.b1, backfill.start.combi.b2

Same as backfill.start.[...].a (where [...] is mono1, mono2, or combi) but for the second potentially simulated trial (suffix .b) of the respective trial type.

n.studies

Positive integer that specifies the number of studies to be simulated, defaults to 1. Due to the long simulation time, it is recommended to first try lower numbers to obtain an estimation of the run-time for larger numbers of studies. Typically, about 1000 studies are recommended to obtain acceptably accurate simulation results.

seed

Optional positive integer that specified the seed to be used for the simulation. The default is sample.int(.Machine$integer.max, 1). Note that reproducibility can only be obtained when the function is executed on exactly the same computing architecture, using identical software versions (e.g. of the compiler, Stan, and R), and the same input specifications. This is due to the internal use of rstan-package for MCMC sampling, which is only reproducible under these restrictions (refer to the Stan reference manual from rstan-package's homepage for more detail).

chains

Optional positive integer that specifies the number of Markov chains to be used during MCMC sampling, defaults to 4. The parameter is given to the method rstan::sampling() from the package rstan-package.

iter

Optional positive integer that specifies the total number of iterations per chain to be used during MCMC sampling, defaults to 6000. The parameter is given to the method rstan::sampling() from the package rstan-package. Be aware that the first warmup iterations among the iter-many samples are discarded. The function enforces iter to be at least 2000, and further that iter is larger or equal than warmup + 1000 (to ensure that at least 1000 samples from each chain are kept).

warmup

Optional positive integer that specifies the number of iterations per chain that are discarded from the total number of iterations, iter. Defaults to 1000. The parameter is given to the method rstan::sampling() from the package rstan-package. The function enforces warmup to be at least 1000, and further that iter is larger or equal than warmup + 1000 (to ensure that at least 1000 samples from each chain are kept).

adapt_delta

Optional numerical that must be at least 0.6 and smaller than 1, defaults to 0.8. The parameter is given to the argument control of the rstan::sampling() method from the package rstan-package. Translates to target Metropolis acceptance probability during Hamiltonian Monte Carlo, respectively NUTS. Can be used to influence the stepsize Stan uses for leapfrog steps during the NUTS algorithm. See rstan::sampling(), rstan-package for details and further reading.

Note: Larger values than the default can be used to reduce the number of Stan warnings on "divergent transition", but will slow down sampling. The function will not permit to set adapt_delta to values below 0.6 (the rstan-package default).

max_treedepth

Optional integer that must be at least 10, defaults to 15. The parameter is given to the argument control of the rstan::sampling() method from the package rstan-package. This is a parameter of the NUTS algorithm. Roughly speaking, NUTS constructs a search tree when generating a new proposal, until a stopping criterion (the NUTS criterion) is satisfied or, alternatively, until the maximum depth of the search tree is reached (to avoid endless looping). The argument max_treedepth allows to control the latter. The maximum treedepth will also automatically constrain the maximum number of leapfrog steps, and should therefore not be below 10. See also rstan::sampling() for details on max_treedepth.

Note: The default maximum treedepth of 15 should usually not be reached during sampling (and never was in tests, although this cannot be excluded in general), otherwise Stan will print a warning. In this case, it is recommended to increase its value.

refresh

Optional integer, defaults to 0. The parameter is given to the rstan::sampling() method from the package rstan-package and controls the amount of intermediate output during calls to rstan::sampling(). Only affects non-parallelized simulations (otherwise it only affects the output of each worker, which is discarded anyway). If a positive value is given, the number is interpreted as the number of samples after which intermediate progress is reported during MCMC sampling (e.g., when set to 1000, Stan prints a progress report every 1000 samples). If refresh is 0 (the default) or negative, no output will be printed during sampling. Changing the value of refresh does not affect the results of the simulation, although it is recommended to not change it to a positive number: the function can lead to many calls to rstan::sampling() (multiple thousand when 1000 trials are simulated), so that intermediate output would not be helpful or interpretable in most settings anyway.

file.name

Optional character string that provides a name for potential output files. No output will be written if file.name or path are not specified, and in this case, all outputs are only returned to R. If both are given, an excel sheet with the simulation results is written to disk.

path

Optional character string that specifies the path to the output directory. No output will be written if file.name or path are not specified, and in this case, all outputs are only returned to R. If both are given and path points to a directory, an excel sheet an excel sheet that contains the simulation results is written to the specified location.

monitor.path

Optional character string that specifies a path to an additional output directory for monitoring simulation progress. More precisely, if monitor.path is not NULL (the default) and a valid file.name was given, the function will write an output file for each completed trial to monitor.path (provided that it specifies a writable directory). This can be used to monitor the progress of (parallelized) simulations, which may take multiple hours to simulate.

working.path

Optional character string that specifies a path to a directory for temporary results. This directory is then used to save MCCM results for a given data scenario under unique (hashed) file names, which allows re-loading MCMC runs whenever the same data scenario occurs in some later trial. This also works for parallelized simulations. File names of temporary results will use the given file.name concatenated with "_tmp" as prefix. Do not modify such files manually to avoid errors. To avoid reloading results from previous simulations (with potentially different MCMC settings or priors), the function will check for files whose name contains the file.name and "_tmp" from the directory specified in working.path before the simulation. If such files are detected, an error will be thrown. To allow the function to automatically delete all such files from the working.path before and after the simulation, consider using the argument clean.working.path.

clean.working.path

Optional logical, defaults to FALSE. Indicates whether the function is allowed to delete any existing files that contain file.name and "_tmp" in their name from the directory specified in working.path. If TRUE, such temporary files will be automatically removed from the working directory before and after simulation. If FALSE, the function will check for such file names initially and throw an error when such files are found.

output.sim.config

Optional logical that specifies whether the input parameters of the call to sim_jointBLRM() are added to the output list (or, if activated, the output files written to disk). Defaults to TRUE. This can be used to save e.g. the seed, dose-toxicity scenarios, and priors of the simulation to allow double-checking the configuration of the simulation at a later point in time.

output.sim.detail

Optional logical, defaults to FALSE. Specifies if a list of more detailed simulation results per simulated study should be created. If so, this list will be written to a separate excel sheet, provided a path is given for outputs. The list contains the number of patients (total, at over/under/target) doses per trial, the number of DLTs per trial, the number of patients per dose for each trial, the number of DLTs per dose for each trial, and the name of the dose selected as MTD for each trial (or the empty string, if the trial did not find an MTD).

Value

List that contains a number of metrics that summarize the results of the simulation for each simulated trial, and, depending on the specification, additional list entries that save the given input specification.

By default, the following list entries are generated for each simulated trial (where [...]

is the corresponding suffix of the trial):

  • ...$'results [...]'
    Overview of number of MTDs per dosing interval (under, target, over), and number of trials that stopped without finding an MTD, either due to the stopping rule (when EWOC-based escalation is used) or due to reaching the maximal patient number.

  • ...$'summary [...]'
    Summary statistics (mean, median, min, max, 2.5% and 97.5% quantile) of the number of patients, patients per dosing interval, DLTs, and DLTs per dosing interval.

  • ...$'MTDs [...]'
    Number of MTDs per dose level and their assumed DLT rates.

  • ...$'#pat [...]'
    Mean, median, minimum and maximum number of patients enrolled at each dose level.

Additionally, if output.sim.config is active, the list will include the following entries:

  • ...$'historical.data'
    Only included when historical data was included in the simulations. Contains the specified cohorts from historical data.

  • ...$'prior'
    The prior distribution for the hyper-parameters that was used by the BLRM.

  • ...$'specifications'
    All other simulation specifications. Includes e. g. reference doses, decision rules, escalation steps, and the seed.

  • ...$'Stan options'
    Options given to Stan, such as the number of MCMC iterations and chains.

Details

The joint BLRM is defined according to (Neuenschwander et al., 2014 and 2016). It allows to perform Bayesian logistic regression to estimate the dose-toxicity relationship of two different monotherapies and combination therapy with these compounds in a joint model, which includes a hierarchical prior for robust borrowing across trials. Refer to the documentation of scenario_jointBLRM() for a detailed model description.

References

Stan Development Team (2020). RStan: the R interface to Stan. R package version 2.21.2. https://mc-stan.org.

Neuenschwander, B., Branson, M., & Gsponer, T. (2008). Critical aspects of the Bayesian approach to phase I cancer trials. Statistics in medicine, 27(13), 2420-2439, doi:10.1002/sim.3230.

Neuenschwander, B., Matano, A., Tang, Z., Roychoudhury, S., Wandel, S., & Bailey, S. (2014). A Bayesian Industry Approach to Phase I Combination Trials in Oncology. In: Zhao. W & Yang, H. (editors). Statistical methods in drug combination studies. Chapman and Hall/CRC, 95-135, doi:10.1201/b17965.

Neuenschwander, B., Roychoudhury, S., & Schmidli, H. (2016). On the use of co-data in clinical trials. Statistics in Biopharmaceutical Research, 8(3), 345-354, doi:10.1080/19466315.2016.1174149.

Babb, J., Rogatko, A., & Zacks, S. (1998). Cancer phase I clinical trials: Efficient dose escalation with overdose control. Statistics in medicine 17(10), 1103-1120.

Zhou, H., Yuan, Y., & Nie, L. (2018). Accuracy, safety, and reliability of novel phase I designs. Clinical Cancer Research, 24(21), 5483-5484 <doi: 10.1158/1078-0432.ccr-18-0168>.

See also

scenario_jointBLRM(), rstan::sampling(), rstan-package.