BOIN : An R Package for Designing Single-Agent and Drug-Combination Dose-Finding Trials Using Bayesian Optimal Interval Designs

This article describes the R package BOIN , which implements a recently developed methodology for designing single-agent and drug-combination dose-ﬁnding clinical trials using Bayesian optimal interval designs (Liu and Yuan 2015; Yuan, Hess, Hilsenbeck, and Gilbert 2016). The BOIN designs are novel “model-assisted” phase I trial designs that can be implemented simply and transparently, similar to the 3 + 3 design, but yield excellent performance comparable to those of more complicated, model-based designs. The BOIN package provides tools for designing, conducting, and analyzing single-agent and drug-combination dose-ﬁnding trials.


Introduction
The objective of a phase I trial is to find the maximum tolerated dose (MTD), which is defined as the dose or dose combination with the dose limiting toxicity (DLT) probability closest to the target DLT rate. Numerous phase I trial designs have been proposed to find the MTD. The 3 + 3 design (Storer 1989) has been dominant in phase I clinical trials for decades due to its simplicity, despite its poor ability to identify the MTD and tendency to treat patients at low doses that are potentially subtherapeutic (Le Tourneau, Lee, and Siu 2009). The 3 + 3 design and its variations are called algorithm-based designs because they use simple, prespecified rules to guide dose escalation. Model-based designs have been proposed to improve the efficiency of identifying the MTD. A typical example of model-based designs is the continual reassessment method (CRM;O'Quigley, Pepe, and Fisher 1990). The CRM yields better performance than the 3 + 3 design, but is statistically and computationally complex, leading practitioners to perceive dose allocations as coming from a "black box", which limits its application in practice.
The Bayesian optimal interval (BOIN) design (Liu and Yuan 2015;Yuan et al. 2016) is a novel phase I trial methodology, belonging to a new class of designs known as model-assisted designs (Yan, Mandrekar, and Yuan 2017;Zhou, Murray, Pan, and Yuan 2018a), that combines the simplicity of algorithm-based designs with the superior performance of model-based designs. The BOIN design can be implemented in a simple way as the 3 + 3 design, but yields excellent performance comparable to that of the more complicated, model-based CRM (Zhou, Yuan, and Nie 2018b;Zhou et al. 2018a). Another example of model-assisted designs is the keyboard design (Yan et al. 2017). Because of their good performance and simplicity, model-assisted designs have been increasingly used in practice. Clertant and O'Quigley (2017) proposed a semiparametric approach that provides the link between the BOIN and other designs, such as the cumulative cohort design (Ivanova, Flournoy, and Chung 2007).
The objective of this article is to provide practitioners a comprehensive, well-documented, easy-to-use R (R Core Team 2020) package to design phase I clinical trials using the BOIN design. The BOIN package (Yuan and Liu 2020) incorporates a number of practical considerations and functionalities that are not considered in the original methodological paper by Liu and Yuan (2015). For example, the BOIN R package allows users to early stop the trial when the number of patients accumulated in a dose reaches a certain cutoff. To the best of our knowledge, this is the first R package that provides a comprehensive implementation of the BOIN design, including simulation functions for generating operating characteristics of the design and implementation functions for conducting a real trial in practice. The existing dose-finding R packages, such as TEQR (Blanchard 2016) and DoseFinding (Bornkamp, Pinheiro, and Bretz 2009;Bornkamp 2019), cover algorithm-based designs (i.e., the cumulative cohort design) and model-based designs (i.e., MCP-Mod design), respectively, and focus on single-agent dose finding. In contrast, the BOIN R package focuses on the model-assisted design, and covers both single-agent and drug-combination dose findings. Drug combination trials are increasingly important in drug development.

Bayesian optimal interval designs 2.1. BOIN design for single-agent trials
Let d 1 < . . . < d J denote the J prespecified doses of the new drug that is under investigation in the trial, p j denote the dose-limiting toxicity (DLT) probability that corresponds to d j , and φ denote the target DLT rate for the MTD. We use n j to denote the number of patients who have been assigned to d j , and y j of them experienced DLTs, j = 1, . . . , J. The BOIN design is illustrated in Figure 1 and summarized as follows: (a) Patients in the first cohort are treated at the lowest or a prespecified dose level.  (b) Letp j = y j /n j be the observed DLT rate at the current dose level j, to assign a dose to the next cohort of patients, • ifp j ≤ λ e , escalate the dose level to j + 1, • ifp j ≥ λ d , de-escalate the dose level to j − 1, • otherwise, i.e., λ e <p j < λ d , stay at the same dose level, where λ e and λ d are prespecified dose escalation and de-escalation boundaries (see Table 1). The method of specifying λ e and λ d is described below.
(c) Repeat step (b) until the prespecified maximum sample size N max is exhausted or the trial is terminated because of excessive toxicity, as described next. At that point, the MTD is selected as the dose for which the isotonic estimate of the DLT probability is closest to the target DLT rate φ. The isotonic estimates of the DLT probabilities can be obtained using the pooled adjacent violators algorithm (Barlow, Bartholomew, Bremner, and Brunk 1972). In the case that the trial is terminated because of excessive toxicity, no dose should be selected as the MTD. For patient safety, during the trial conduct, the BOIN design imposes a dose elimination/safety stopping rule as follows: Suppose j is the current dose level. If P(p j > φ|y j , n j ) > 0.95 and n j ≥ 3, dose level j and higher doses are eliminated from the trial. The trial is terminated if the lowest dose is eliminated.
The posterior probability P(p j > φ|y j , n j ) can be evaluated on the basis of a beta-binomial model, assuming y j | n j , p j ∼ Binom(n j , p j ) and p j ∼ Beta(1, 1) ≡ Unif(0, 1). Under the beta-binomial model, the posterior distribution of p j arises as The dose escalation and de-escalation boundaries (λ e , λ d ) are obtained by minimizing the chance of making incorrect dose escalation and de-escalation (Liu and Yuan 2015). To do that, the BOIN design requires the investigator(s) to specify φ 1 , which is the highest DLT rate that is deemed to be underdosing such that dose escalation is required, and φ 2 , which is the lowest DLT rate that is deemed to be overdosing such that dose de-escalation is required. Liu and Yuan (2015) provided general guidance to specify φ 1 and φ 2 , and recommended default values of φ 1 = 0.6φ and φ 2 = 1.4φ for general use. When needed, the values of φ 1 and φ 2 can be calibrated to achieve a particular requirement of the trial at hand. For example, if more conservative dose escalation is required, setting φ 2 = 1.2φ may be appropriate. We should avoid setting the values of φ 1 and φ 2 very close to φ. This is because the small sample sizes of typical phase I trials prevent us from differentiating the target DLT rate from the probabilities that are close to it. For example, at the significance level of 0.1, there is only 7% power to distinguish 0.25 from 0.35 with 15 patients at each dose, based on Fisher's exact test. Given φ 1 and φ 2 , the optimal escalation and de-escalation boundaries (λ e , λ d ) that minimize the incorrect decision of dose escalation and de-escalation arise as . Table 1 provides the dose escalation and de-escalation boundaries (λ e , λ d ) for commonly used target DLT rates φ using the recommended default values φ 1 = 0.6φ and φ 2 = 1.4φ.
For example, given the target DLT rate φ = 0.25, the corresponding escalation boundary λ e = 0.197 and the de-escalation boundary λ d = 0.298, that is, escalate the dose if the observed DLT rate at the current dosep j ≤ 0.197 and de-escalate the dose ifp j ≥ 0.298. Liu and Yuan (2015) showed that (λ e , λ d ) are the boundaries corresponding to the likelihood ratio test and Bayes factor, and thus the resulting BOIN design is optimal with desirable finite-sample and large-sample properties, i.e., long-memory coherence and consistency.
As the observed DLT ratep j is the most natural and intuitive estimate of p j that is accessible by non-statisticians, it is easy to explain the BOIN design to clinicians. In addition, thanks to the feature that the BOIN design guarantees de-escalating the dose whenp j is higher than the de-escalation boundary λ d , it is particularly easy for clinicians and regulatory agents to assess the safety of a trial using the BOIN design. For example, given a target DLT rate φ = 0.25, we know a priori that a phase I trial using the BOIN design guarantees de-escalating the dose if the observed DLT rate is higher than 0.298 (with the default values of φ 1 and φ 2 ). Accordingly, the BOIN design also allows users to easily calibrate the design to satisfy a specific safety requirement mandated by regulatory agents through choosing an appropriate target DLT rate φ or φ 2 . For example, suppose for a phase I trial with a new compound, the regulatory agent mandates that if the observed DLT rate is higher than 0.25, the dose must be de-escalated. We can easily fulfill that requirement by setting the target DLT rate φ = 0.21, under which the BOIN automatically guarantees de-escalating the dose if the observed DLT ratep j > λ d = 0.250. If needed, the de-escalation boundary λ d can be further fine tuned by calibrating the value of φ 2 . Such flexibility and transparency gives the BOIN design an important advantage over many existing designs in practice.

BOIN designs for drug-combination trials
Drug combination therapy provides a critical approach to improve treatment efficacy and overcome resistance to monotherapies. One challenge in designing drug-combination trials is that, due to complicated drug-drug interactions, combinations are only partially ordered in toxicity probabilities. Consider a trial combining J doses of drug A, denoted as A 1 < A 2 < . . . < A J , and K doses of drug B, denoted as B 1 < B 2 < . . . < B K . Let A j B k denote the combination of A j and B k , and p jk denote the true DLT rate of A j B k . We know a priori that A j B k is more toxic than A j B k for j > j and k > k, but typically do not know their toxicity order if j > j but k < k. Another challenge in designing combination trials is the existence of the MTD contour in the two-dimensional dose space, as shown in Figure 2. In other words, multiple MTDs may exist in the J × K dose matrix. Depending on the application, combination trials may be interested in finding a single MTD or the MTD contour (i.e., multiple MTDs; Yuan and Zhang 2017). Because of these challenges, the singleagent phase I trial designs cannot be directly used for combination trials. In what follows, we first describe the BOIN drug-combination design (Lin and Yin 2017a) that aims to find a single MTD in the two-dimensional dose combination matrix, followed by a BOIN-waterfall design (Zhang and Yuan 2016) that aims to find the MTD contour (i.e., multiple MTDs)

BOIN drug-combination design for finding a single MTD
Letp jk = y jk /n jk denote the observed DLT rate at dose combination A j B k , where n jk denotes the number of patients treated at A j B k and y jk denotes the number of patients who experienced DLT at A j B k . Given the current dose combination A j B k , define an admissible dose escalation set as A E = {A j+1 B k , A j B k+1 } and an admissible dose de-escalation set as The BOIN drug-combination design makes the decision of dose escalation/de-escalation based on the same rule as the single-agent BOIN design described previously. That is, escalate the dose ifp jk ≤ λ e , and de-escalate the dose ifp jk ≥ λ d . The only difference is that, in combination trials, when we decide to escalate or de-escalate the dose, there is more than one neighbor dose to which we can move. For example, when we escalate/de-escalate the dose, we can escalate/de-escalate either the dose of drug A or the dose of drug B. The BOIN drug-combination design makes this choice based on P(p jk ∈ (λ e , λ d )|D jk ), which measures how likely a dose combination is located within the acceptable toxicity interval (λ e , λ d ), where D jk = (y jk , n jk ). Specifically, the BOIN drug-combination design can be described as follows: (a) Patients in the first cohort are treated at the lowest dose combination A 1 B 1 or a prespecified dose combination.
(b) Suppose the current cohort is treated at dose combination A j B k . Proceed as follows to assign a dose to the next cohort of patients.
• Ifp jk ≤ λ e , we escalate the dose to the combination that belongs to A E and has the largest value of P{p j k ∈ (λ e , λ d )|D jk }.
• Ifp jk ≥ λ d , we de-escalate the dose to the combination that belongs to A D and has the largest value of P{p j k ∈ (λ e , λ d )|D jk }.
• Otherwise, i.e., λ e <p jk < λ d , the dose stays at the same combination A j B k . (c) Step (b) is continued until the maximum sample size is reached or the trial is terminated because of excessive toxicity.
During dose escalation and de-escalation, if the two combinations in A E or A D have the same value of P{p j k ∈ (λ e , λ d )|data}, we randomly choose one with equal probability. If no dose combination exists in the sets of A E and A D (i.e., we are at the boundaries of the dose matrix), we retain the current dose combination. After the trial is completed, the MTD is selected as the dose combination with the estimated DLT rate closest to φ. The estimates of DLT rates are obtained using isotonic regression as described previously, but in a matrix form. More details on the BOIN drug-combination design can be found in Lin and Yin (2017a).

BOIN-waterfall design for finding the MTD contour
Because of the existence of the MTD contour and the fact that doses on the MTD contour may have different efficacy due to drug-drug interactions, for many drug combination trials, it is of intrinsic interest to find multiple MTDs (or the MTD contour). The efficacy of the MTDs can be evaluated in subsequent phase II trials or simultaneously in phase I-II trials. Zhang and Yuan (2016) extended the BOIN design to achieve this design goal. For convenience, we refer to that design as the BOIN-waterfall design. As illustrated in Figure 3, the BOIN-waterfall design takes the divide-and-conquer strategy to divide the task of finding the MTD contour into a sequence of one-dimensional dose-finding tasks, known as "subtrials" (Yuan and Yin 2008). These subtrials are conducted sequentially from the top of the matrix to the bottom.
Specifically, the BOIN-waterfall design can be described as follows: 1. Divide the J × K dose matrix into J subtrials S J , . . . , S 1 , according to the dose level of drug A: . . .
Note that subtrial S J also includes lead-in doses A 1 B 1 , A 2 B 1 , . . . , A J B 1 (the first column of the dose matrix) to impose the practical consideration that the trial starts at the lowest dose combination. Within each subtrial, the doses are fully ordered with monotonically increasing toxicity.
2. Conduct the subtrials sequentially using the single-agent BOIN design as follows: (i) Conduct subtrial S J , starting from the lowest dose combination A 1 B 1 , to find the MTD. We call the dose selected by the subtrial the "candidate MTD" to highlight that the dose selected by the individual subtrial may not be the "final" MTD that we will select at the end of the trial. The final MTD selection will be based on the data collected from all the subtrials. The objective of finding the candidate MTD is to determine which subtrial will be conducted next and the corresponding starting dose combination.
(ii) Assuming that the current subtrial S j , j = J, . . . , 2, selects dose A j * B k * as the candidate MTD, next, conduct subtrial S j * −1 with the starting dose A j * −1 B k * +1 . That is, the next subtrial to be conducted is the one with the dose of drug A that is one level lower than the candidate MTD found in the previous subtrial.
After identifying the candidate MTD of subtrial S j * −1 , the same rule is used to determine the next subtrial and its starting dose combination. See Figure 3 for an example. Figure 3: Illustration of the waterfall design for a combination trial of drug A with 3 levels and drug B with 4 levels. The doses in the rectangle form a subtrial, and the asterisk denotes the candidate MTD. As shown in panel (a), the trial started by conducting the first subtrial with the starting dose A 1 B 1 . After the first subtrial identified A 3 B 2 as the candidate MTD, we conducted the second subtrial with the starting dose A 2 B 3 (see panel (b)). After the second subtrial identified A 2 B 3 as the candidate MTD, we conducted the third subtrial with the starting dose A 1 B 4 (see panel (c)). All subtrials were conducted using the single-agent BOIN design. After all the subtrials were completed, we selected the MTD contour based on the data from all the subtrials, as shown in panel (d).
3. Estimate the DLT rate p jk based on the toxicity data collected from all the subtrials using matrix isotonic regression. For each row of the dose matrix, select the MTD as the dose combination that has the estimate of DLT rate closest to the target DLT rate φ, unless all combinations in that row are overly toxic.
For the BOIN-waterfall design, one key issue is to determine when to end the current subtrial and initiate the next one. Zhang and Yuan (2016) recommended the following stopping rule to switch from one subtrial to the next one: At any time during the subtrial, if the total number of patients treated at the current dose reaches a certain prespecified number of patients, say n stop , we stop the subtrial, select the candidate MTD and initiate the next subtrial.
The rationale for the stopping rule is that when patient allocation concentrates at a given dose, it indicates that the dose finding might have converged to the MTD, and thus we can stop the trial and claim the MTD. This stopping rule allows the sample size of the subtrials to be automatically adjusted according to the difficulty of the dose finding. Another attractive feature of the above approach is that it automatically ensures that a certain number of patients are treated at the MTD. Conventionally, we achieve this by adding cohort expansion after identifying the MTD. In practice, we recommend n stop > 9 to ensure reasonable operating characteristics. Although the above stopping rule provides an automatic, reasonable way to determine the sample size for a subtrial, in some cases, it is desirable to put a cap on the maximum sample size of the subtrials. This can be done by adding an extra stopping rule as follows: Stop the subtrial if its sample size reaches N max , where N max is the prespecified maximum sample size.

The R package BOIN
The R package BOIN is freely available from the Comprehensive R Archive Network (CRAN) at https://CRAN.R-project.org/package=BOIN. It contains functions for implementing the BOIN design for both single-agent and drug-combination dose-finding trials (see Table 2). We proceed with a small tutorial to present the use of the package. As an illustration, we load the package: R> install.packages("BOIN") R> library("BOIN")

Single-agent trials
To design and conduct a single-agent dose-finding trial, we follow 3 steps as described below, which involve functions get.boundary(), get.oc() and selet.mtd(), respectively. The objects (or results) returned by these functions can be displayed in a user friendly form by using the functions summary(), print(), and plot().

. Obtain dose escalation and de-escalation boundaries
Like the algorithm-based design, the dose escalation and de-escalation rule of the BOIN design can be determined before the onset of the trial, making its implementation transparent and simple. To obtain the dose escalation and de-escalation boundaries (λ e , λ d ), we run the function: get.boundary(target, ncohort, cohortsize, n.earlystop = 100, p.saf = 0.6 * target, p.tox = 1.4 * target, cutoff.eli = 0.95, extrasafe = FALSE, offset = 0.05) This function has the following arguments: • target: The target DLT rate.
• ncohort: The total number of cohorts.

Function Description
Object returned get.boundary Generate dose escalation and deescalation boundaries.
A list including dose escalation and de-escalation boundaries.
select.mtd Select the maximum tolerated dose (MTD) for single-agent trials.
A list including the selected MTD and estimated toxicity probability for each dose.
get.oc Generate operating characteristics for single-agent trials.
A list including the selection percentage, the number of patients treated at each dose, and average total number of patients.

next.comb
Determine the dose combination for treating the next cohort of new patients for drug-combination trials that aim to find an MTD.
A list including the recommended dose combination for the next cohort of patients.
next.subtrial Determine the starting dose combination and the dose-searching space for the next subtrial in BOIN-waterfall design.
A list including dose combinations included in the subtrial and the starting dose combination.
select.mtd.comb Select the MTD or MTD contour for drug combination trials.
A list containing the selected MTD or MTD contour, and estimated toxicity probability for each dose combination.
get.oc.comb Generate operating characteristics for drug-combination trials.
A list containing the selection percentage, the number of patients treated at each dose, and average total number of patients.
summary Generate summary for the objects returned by other functions.
Descriptive results printed to the console.
plot Generate flowchart of the design and barplot for simulation results.

Figures.
print Print objects returned by other functions.
Objects returned from other functions. Table 2: Overview of the user visible functions in BOIN. Please consult the documentation (e.g., help("get.boundary")) for function arguments and detailed return types.
• n.earlystop: The early stopping parameter. If the number of patients treated at the current dose reaches n.earlystop, stop the trial early and select the MTD based on the observed data. The default value of n.earlystop = 100 essentially turns off this type of early stopping.
• p.saf: The highest toxicity probability that is deemed subtherapeutic (i.e., below the MTD) such that dose escalation should be made. The default value of p.saf = 0.6 * target.
• p.tox: The lowest toxicity probability that is deemed overly toxic such that dose deescalation is required. The default value of p.tox = 1.4 * target.
• cutoff.eli: The cutoff to eliminate the overly toxic dose for safety. We recommend the default value cutoff.eli = 0.95 for general use.
• offset: A small positive number (between 0 and 0.5) to control how strict the stopping rule is when extrasafe = TRUE. A larger value leads to a stricter stopping rule. The default value offset = 0.05 generally works well.

R> plot(bound)
A few remarks are in order. First, the flowchart generated by plot() provides a transparent description of the trial design that can be easily understood by clinicians. We strongly recommend including it in the trial protocol. Second, as described previously, in practice, we should avoid setting the values of p.saf and p.tox very close to the target because the small sample sizes of typical phase I trials prevent us from being able to discriminate the target DLT rate from the rates close to it. The default values φ 1 = 0.6φ and φ 2 = 1.4φ are strongly recommended and generally yield superior operating characteristics.
Lastly, the BOIN design has two built-in stopping rules: (i) stop the trial if the lowest dose is eliminated due to toxicity. In this case, no dose should be selected as the MTD; (ii) stop the trial and select the MTD if the number of patients treated at the current dose reaches n.earlystop. The first stopping rule is a safety rule to protect patients from the case in which all doses are overly toxic. The rationale for the second stopping rule is that when the number of patients assigned to a dose is large (i.e., reaches n.earlystop), this means that the dosefinding algorithm has approximately converged. Thus, we can stop the trial early and select the MTD to save sample size and reduce the trial duration. The default value n.earlystop = 100 essentially turns off this type of early stopping rule. Setting n.earlystop at a value like 12 can potentially save sample size and finish the trial early. The trade-off is that it may affect the MTD selection percentage and decrease the rate of stopping for safety if the first dose is overly toxic. The value of n.earlystop should be calibrated by simulation to obtain desirable operating characteristics. In general, we recommend n.earlystop = 9 to 18. Our experience is that this stopping rule is particularly useful when there is strong prior knowledge that the first dose is safe, since a major side effect of using the stopping rule is that it decreases the rate of stopping for safety when the first dose is actually overly toxic.
Although the BOIN design has a built-in safety stopping rule (i.e., stopping rule (i) described above), for some applications, investigators may prefer a stricter stopping rule for extra safety when the lowest dose is possibly overly toxic. Setting extrasafe = TRUE imposes the following stronger stopping rule: Stop the trial if (1) the number of patients treated at the lowest dose ≥ 3, and (2) P(DLT rate of the lowest dose > target | data) > cutoff.eli -offset.
As a trade-off, the stricter stopping rule will decrease the MTD selection percentage when the lowest dose actually is the true MTD. When using the option extrasafe = TRUE, we recommend the default value offset = 0.05, but users can calibrate the value of offset to obtain desired operating characteristics. In practice, offset is rarely greater than 0.2. If we set extrasafe = TRUE to turn on the extrasafe feature, the output will include the extra stopping boundaries, as follows, R> bound2 <-get.boundary(target = 0.3, ncohort = 10, cohortsize = 3, + extrasafe = TRUE) R> summary(bound2) ......
In addition to the default stopping rule (i.e., stop the trial if the lowest dose is eliminated), the following more strict stopping safety rule will be used for extra safety: stop the trial if (1) the number of patients treated at the lowest dose >= 3, AND (2) Pr(the DLT rate of the lowest dose > 0.3 | data) > 0.9, which corresponds to the following stopping boundaries: The Step 2. Simulate operating characteristics For protocol preparation, it is often useful to obtain the operating characteristics of the design. The function get.oc() can be used for this purpose.
• startdose: The starting dose level for treating the first cohort of patients. The default value is startdose = 1, i.e., starting from the lowest dose.
• titration: If titration = TRUE, dose titration is performed to accelerate dose escalation at the beginning of the trial, where patients are treated one by one (i.e., cohort size = 1), starting from startdose. If no DLT is observed, escalate the dose; otherwise switch to the specified cohort size = cohortsize.
• ntrial: The number of trials to be simulated.
• seed: Set a seed for the random number generator.
Dose titration, by setting titration = TRUE, accelerates dose escalation and improves the selection percentage of the MTD when the number of patients is small (e.g., ncohort/(number of doses) < 1.5) and low doses are expected to be safe. The tradeoff is that if the low doses are toxic, the dose titration increases the risk of overdosing patients.
Using the same setting as above and assuming that the true toxicity scenario is p.true = (0.05, 0.15, 0.30, 0.45, 0.6), here we show how to obtain the operating characteristics based on 1000 simulated trials.

R> plot(oc.single)
Step 3. Select the MTD When the trial is completed, based on the observed data, we can select the MTD using the function select.mtd().
• ntox: The vector recording the number of patients who experienced toxicity at each dose level. Arguments cutoff.eli, extrasafe and offset are the same as (and should be consistent with) those in functions get.boundary() and get.oc(). Assume that when the trial is completed, the number of patients treated at the five doses is n = (3, 3, 15, 9, 0) and the corresponding number of patients who experienced toxicity is y = (0, 0, 4, 4, 0).
We can plot the estimates of DLT rates ( Figure 5) by calling the plot() function.

Drug-combination trials aiming to find a single MTD
Designing a drug-combination trial follows the similar three steps as those for single-agent trials. This involves four functions get.boundary(), next.comb(), get.oc.comb() and selet.mtd.comb().

Step 1. Determine dose escalation and de-escalation
The BOIN drug-combination design uses the same dose escalation and de-escalation boundaries as the single-agent designs, which can be obtained using get.boundary() as described previously. For drug combination trials, the challenge is that when the number of observed DLTs crosses the dose escalation and de-escalation boundaries, we need to determine which drug (i.e., drug A or drug B) to be escalated or de-escalated. The function next.comb() is used to make such a decision.
next.comb(target, npts, ntox, dose.curr, n.earlystop = 100, p.saf = 0.6 * target, p.tox = 1.4 * target, cutoff.eli = 0.95, extrasafe = FALSE, offset = 0.05) This function takes the data from patients who have been enrolled into the trial as the input, and outputs the dose combination for treating the next cohort of new patients. It shares a similar set of arguments with the function get.boundary() described previously, with three additional arguments: • npts: The matrix recording the number of patients treated at each dose combination.
• ntox: The matrix recording the number of patients who experienced toxicity at each dose combination.
• dose.curr: The current dose combination, i.e., the dose combination that was used to treat the most recently enrolled cohort of patients.
Therefore, we should de-escalate the dose and treat the third cohort of patients at dose (1, 1). We repeat this procedure until the maximum sample size is reached.
Step 2. Simulate operating characteristics The function get.oc.comb() can be used to obtain the operating characteristics of the BOIN drug-combination design. This function shares the same set of arguments as the function get.oc() described previously (for single-agent trials), except that p.true is now a matrix (rather than a vector), startdose is a vector of length 2 (rather than a scalar), and an additional arguments mtd.contour is used to indicate whether we are interested in finding the MTD or the MTD contour. To find a single MTD, we should set mtd.contour = FALSE. Again, performing dose titration, by setting titration = TRUE, is useful to accelerate dose escalation and improve the selection percentage of the MTD when the number of patients is small (e.g., ncohort/(number of doses) < 1.5) and low doses are expected to be safe. However, when the low doses are toxic, dose titration may increase the risk of overdosing patients. We can use the plot() function to plot the figure illustrating dose escalation/de-escalation and operating characteristics such as selection percentage (Figure 6), average number of patient treated, and average number of toxicity observed at each dose combination.

Step 3. Select an MTD
When the trial is completed, based on the observed data, we can select the MTD using the function select.mtd.comb().

Drug-combination trials aiming to find the MTD contour
Designing a drug-combination trial for finding the MTD contour also follows the similar three steps as those for single-agent trials. This involves four functions get.boundary(), next.subtrial(), get.oc.comb() and selet.mtd.comb().

. Obtain dose escalation and de-escalation boundaries and determine subtrials
The BOIN-waterfall design is used to find the MTD contour, which divides the task of finding the MTD contour into a sequence of subtrials. Each subtrial is conducted using the BOIN single-agent design, whose dose escalation and de-escalation boundaries can be obtained using get.boundary() as described previously. The function next.subtrial() is used to determine the next subtrial based on the observed data.
next.subtrial(target, npts, ntox, p.saf = 0.6 * target, p.tox = 1.4 * target, cutoff.eli = 0.95, extrasafe = FALSE, offset = 0.05) As described in Figure 3 and Zhang and Yuan (2016), the first subtrial is prespecified and starts from the lowest dose (1, 1). When the first subtrial is completed, we call the function next.subtrial() to determine the starting dose combination and the dose-searching space for next subtrial. The function next.subtrial() shares a similar set of arguments with the function get.boundary() described previously, with two additional arguments: • npts: The matrix recording the number of patients treated at each dose combination.
• ntox: The matrix recording the number of patients who experienced toxicity at each dose combination.
Step 2. Simulate operating characteristics The function get.oc.comb() can be used to obtain the operating characteristics of the BOINwaterfall design. This function shares the same set of arguments as the function get.oc() described in Section 3.1 (for single-agent trials), except that p.true is now a matrix (rather than a vector), startdose is a vector of length 2 (rather than a scalar), and an additional argument mtd.contour is used to indicate whether we are interested in finding an MTD or the MTD contour. R> p.true <-matrix(c(0.01, 0.03, 0.10, 0.20, 0.30, 0.03, 0.05, 0.15, 0.30, + 0.60, 0.08, 0.10, 0.30, 0.60, 0.75), nrow = 3, ncol = 5, byrow = TRUE) R> oc.comb2 <-get.oc.comb(target = 0.3, p.true, ncohort = c(10, 5, 5)  We can use the plot() function to plot the figure illustrating the waterfall design (similar to Figure 3) and operating characteristics such as selection percentage (Figure 7), average number of patient treated, and average number of toxicity observed at each dose combination.
R> plot(oc.comb2) Step 3. Select the MTD contour When the trial is completed, based on the observed data, we can select the MTD contour using the function select.mtd.comb(), described previously. Since we are interested in finding the MTD contour, we should set argument mtd.contour = TRUE.

Single-agent phase I trial
Consider a single-agent phase I trial with 5 dose levels, in which the objective is to find the MTD with a target DLT rate of 0.3. The maximum sample size is 30 patients, treated in cohort sizes of 3. To design and conduct this trial, we first ran the following commands: R> bound <-get.boundary(target = 0.3, ncohort = 10, cohortsize = 3) R> summary(bound) R> plot(bound) This yields the dose escalation and de-escalation boundaries as shown in Table 3, and a flowchart of the trial design similar to Figure 1. The trial started by treating the first cohort of 3 patients at dose level 1 and none of the patients had dose limiting toxicity (DLT). According to the dose escalation and de-escalation rule provided in Table 3, we escalated the dose to level 2 to treat the second cohort of 3 patients, none of whom experienced DLT. Thus, we escalated the dose to level 3 and treated the third cohort of patients, 2 of whom experienced DLT. Based on Table 3, we de-escalated the dose back to level 2 and treated the fourth cohort of patients, one of whom experienced DLT. We then escalated the dose to level 3 and treated the fifth cohort of patients, none of whom experienced DLT. Therefore, q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q the sixth cohort was also treated at dose level 3. Figure 8 shows the dose assignment for all 30 patients. At the end of the trial, the number of patients and the number of DLTs at the 5 doses were n = c(3, 6, 18, 3, 0) and y = c(0, 1, 5, 3, 0). To select the MTD, we ran the following commands: R> n <-c(3, 6, 18, 3, 0) R> y <-c(0, 1, 5, 3, 0) R> sel.single <-select.mtd(target = 0.3, ntox = y, npts = n) R> summary(sel.single) The MTD is dose level 3 ----(----,----) ----NOTE: no estimate is provided for the doses at which no patient was treated.
which recommended the dose level 3 as the MTD, with the estimated DLT rate = 0.28 and the 95% CI = (0.10, 0.50). And the illustration of the single-agent phase I trial is shown in Figure 8.

Drug-combination trial to find a single MTD
Consider a drug-combination trial that combines 3 doses of drug A and 5 doses of drug B. The objective is to find an MTD with a target DLT rate of 0.3. The maximum sample size is 30 patients, treated in cohort size of 3. The trial started by treating the first cohort of 3 patients at the lowest dose combination (1, 1), at which no DLT was observed. The observed data were n = where n records the number of patients treated at each dose combination, and y records the number of patients who experienced DLT at each dose combination. In matrices y and n, entry (j, k) records the data associated with dose combination (j, k). To determine the dose for the second cohort of patients, we used the following code: R> set.seed(4) R> n <-matrix(c(3, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0), nrow = 3, + ncol = 5, byrow = TRUE) R> y <-matrix(c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0), nrow = 3, + ncol = 5, byrow = TRUE) R> nxt.comb <-next.comb(target = 0.3, npts = n, ntox = y, + dose.curr = c(1, 1)) R> summary(nxt.comb) The recommended dose combination for the next cohort of patients is (1, 2). which recommended to escalate the dose to combination (1, 2). Therefore, we treated the second cohort of patients at dose combination (1, 2). In the second cohort, no patient experienced DLT, so the updated data matrices became To determine the dose for the third cohort of patients, we again called next.comb() with updated y, n and dose.curr.
The function recommended escalating the dose to (2, 2) for treating the third cohort of patients. We repeated this procedure until the maximum sample size was reached. Figure 9 shows the dose assignments for all 30 patients. For example, at dose combination (3, 4) when completing the eighth cohort, there were two DLTs, based on the accumulating toxic information on this dose combination level, the function recommended de-escalating the dose to combination (3, 3). When the trial was completed, the number of patients treated at each dose combination and the corresponding number of patients who experienced toxicity at each dose combination were We called the function select.mtd.comb(), which recommended dose combination (3, 3) as the MTD.

Drug-combination trial to find the MTD contour
Consider a drug-combination trial that combines 3 doses of drug A and 5 doses of drug B. The objective is to find the MTD contour (multiple MTDs) with a target DLT rate of 0.25. The trial started with the first subtrial, which consisted of 7 ordered combinations This subtrial was conducted using the BOIN design in a way similar to those as described in Example 4.1. The starting dose for this subtrial was A 1 B 1 and n.earlystop was set as 12 such that the subtrial stopped when the number of patients treated at that dose reached 12. The first subtrial resulted in the following data: Based on the data, we called the function next.subtrial() to obtain the doses for the next subtrial.

Discussion
As a model-assisted design, the BOIN design combines the simplicity of the algorithm-based design and the superior performance of the model-based design, providing an attractive approach to designing phase I clinical trials. The R package BOIN provides an easy-to-use and well-documented tool to implement the BOIN designs for single-agent and drug-combination phase I trials. The package allows users to generate the dose escalation and de-escalation boundaries, plot the flowchart, and simulate the operating characteristics of the design for trial protocol preparation. It also provides functions to allow users to make the real-time decision of dose assignment and MTD selection for conducting a real trial. By facilitating the use of novel adaptive designs, the BOIN package has potential to improve the efficiency of phase I trials and accelerate drug development. Recently, the BOIN design has been extended to account for toxicity grade, continuous, and semi-continuous toxicity endpoints (Mo, Yuan, Xu, Mandrekar, and Yin 2018), late-onset toxicity (Yuan, Lin, Li, Nie, and Warren 2018), and phase I-II trials (Lin and Yin 2017b;Takeda, Taguri, and Morita 2018). We plan to incorporate these new developments into the future version of the BOIN package.