approach since the first chemotherapies' development in the late 1940s. [2][3][4] There are several types of IT, including monoclonal antibodies, non-specific immunotherapies, oncolytic virus therapy, T-cell therapy, and cancer vaccines. IT has already led to major treatment breakthroughs for several cancers, including brain cancer, breast cancer, bladder cancer, and other types of cancer. 5 The emergence of IT has challenged the assumption that both efficacy and toxicity response rates increase monotonically with the dose, which most of the phase I dose-finding trials rely on to identify the maximum tolerated dose (MTD) solely based on toxicity outcomes. Indeed, IT enhances the immune system's innate power to stop the growth of cancer cells and prevent cancer from spreading to other parts of the body. Consequently, IT is not necessarily administered at its MTD to achieve an optimal therapeutic effect, making the conventional phase I trial design targeting MTD inappropriate for IT trials. Instead, the primary goal of a dose-finding trial for IT is to identify the optimal biological dose (OBD), which achieves the optimal overall therapeutic effect among all the candidate doses. To find the OBD for IT requires monitoring both the toxicity and efficacy outcomes in a single trial, which is typically referred to as the phase I/II clinical trial. 6 Due to different biological mechanisms of cytotoxic agents (e.g., radiotherapy and chemotherapy) and IT, the definitions of the efficacy responses also differ. In general, a patient responds favorably to a cytotoxic agent if he/she has achieved at least partial remission of the tumor after receiving the treatment. On the other hand, because IT aims to stop the growth of tumors rather than kill the cancer cell, stable disease is often considered a positive signal, and only the occurrence of disease progression is treated as a negative efficacy response. 7 As a result, the progression-free survival rate rather than the objective tumor response rate is often used as the efficacy outcome in the IT trial.
Using the progression-free survival rate as the response outcome in a phase I/II dose-finding trial causes many essential and practical issues. First of all, patients developing disease progression in oncology trials are typically seriously ill. Indeed, many patients will develop cancer-related death sooner after experiencing disease progression. Even if not, they will be treated off the trial with other second-line drugs because it is highly unethical to continue treating a seriously ill patient with an ineffective drug. Consequently, if the toxicity event does not happen before the disease progression, it will be censored at the disease progression time. In other words, the happening of disease progression will terminate the happening of toxicity, but not vice versa. In the statistical community, this kind of data is often referred to as the semi-competing risks data. 8 Secondly, the evaluation of progression-free survival rate requires a relatively long follow-up time compared with evaluating the objective tumor response rate, which may cause the late-onset outcome issue in addition to the semi-competing risks. The late-onset happens when at least part of the patients' response outcomes are missing due to incomplete follow-up at the interim analysis time and maybe ascertainable long after that time. The missing data caused by late-onset in a dose-finding trial are not-missing-at-random and should be appropriately handled to avoid biased statistical inference. 9 The naive solution of temporarily suspending the trial to allow all enrolled patients to be fully evaluated before a new dose assignment is often unrealistic because it may result in an unfeasibly long trial, wastes resources, and causes a tremendous administrative burden. Besides, it is ethically undesirable to delay a new patient's treatment while waiting for previous patients' outcomes.
We develop a novel Bayesian adaptive dose-finding design accounting for semi-competing risks outcomes for IT trials in the presence of late-onset issue. We first decompose the joint toxicity-efficacy probability into the product of marginal and conditional probabilities and then re-construct the likelihood function based on each patient's actual followup time. We also develop a data augmentation method to efficiently draw the posterior samples for the parameters of interest from a series of beta-binomial distributions, 10 which substantially reduces the computational burden. We propose a curve-free dose-finding algorithm to adaptively identify the OBD based on a toxicity-efficacy tradeoff function, making no parametric model assumptions on the dose-response relationships.
Our research is motivated by a phase I/II IT trial conducted at Indiana University Melvin and Bren Simon Comprehensive Cancer Center. This trial aims to find the OBD for a novel programmed cell death protein 1 (PD-1) immune checkpoint inhibitor for patients with acute myeloid leukemia (AML) and high-risk myelodysplasia. This inhibitor acts to inhibit the association of the programmed death-ligand 1 with its receptor PD-1, which is involved in suppressing the patient's immune system. 11 Five dose levels (0.1, 0.3, 0.5, 0.8, 1.0 mg/kg) of the inhibitor will be investigated, and the prepared doses will be administered by slow injection over 10 min. A maximum of 60 patients will be accrued to the trial. The primary efficacy outcome is the 1-year progression-free survival rate. The co-primary toxicity outcome is the grade 3 and 4 hematological and non-hematological toxicity rate associated with treatment, which will also be ascertained in a 1-year period. Any patients experiencing disease progression will be treated off the trial with the second-line treatment azacitidine, confounding the toxicity outcome's evaluation. Hence, the issue of semi-competing risks arises. In addition, a new cohort of patients will enter the trial every 6 months for dose allocation before all the toxicity and efficacy outcomes become available for previous patients in the trial. Therefore, this IT trial is also subject to the late-onset issue. Designing a rigorous and efficient dose-finding trial handling late-onset semi-competing risks toxicity-efficacy outcomes is the challenge of this study.
Many phase I/II clinical trial designs for IT have been proposed in the literature. For categorical outcome, Messer et al. 12 proposed a toxicity evaluation design for IT trials, which is based on safety hypothesis testing and an algorithm similar to a 3 + 3 design. Liu and Johnson 13 developed a robust phase I/II dose-finding design that modeled the toxicity and efficacy using a flexible Bayesian dynamic model and borrowed information across doses without imposing stringent parametric assumptions on the shape of the dose-toxicity and dose-efficacy curves. Lin and Yin 14 proposed the STEIN design as the first toxicity-efficacy interval design for optimal dose-finding trials considering both efficacy and toxicity outcomes. Liu et al. 15 developed a latent variable model to capture the correlations of the toxicity and efficacy outcomes in a phase I/II trial. Zhou et al. 16 developed a utility-based Bayesian optimal interval (U-BOIN) design to identify the OBD. Han et al. 17 developed a two-stage nonparametric dose-finding design (TSNP) for IT trial. For lateonset outcomes, Liu et al. 18 treated the late-onset outcomes as missing data and proposed a data augmentation algorithm to impute the missing data under the Bayesian framework. The work was further extended by Jin et al. 19 to the phase I/II trial design considering both toxicity and efficacy outcomes. Zhang and Zang 20 proposed a conditional weighted likelihood (CWL) method to address the jointly late-onset toxicity-efficacy outcomes issue without making any parametric model assumption. To the best of our knowledge, this paper presents the first phase I/II clinical trial design for IT that accounts for semi-competing risks late-onset outcomes.
We have previously developed another Bayesian adaptive design for phase I/II clinical trials with late-onset outcomes. 21 The differences between the previous one and the new design proposed in this paper are listed as follows.
(1) We consider the competing-risk late-onset outcomes in the previous design where either the happening of disease progression or dose-limiting toxicity will terminate the follow-up for another event. In this paper, we focus on the semicompeting risk late-onset outcomes such that only the happening of disease progression can early terminate a patient's follow-up. ( 2) We use a piece-wise cause-specific hazard model to characterize the late-onset competing risk outcome in the previous paper and propose to use the uniform distribution to approximate the late-onset semi-competing risk outcome in this paper. (3) We use a parametric continuation-ratio to capture the correlation between the toxicity and efficacy outcomes in the previous design and use the strategy of decomposing the joint distribution into the marginal and conditional distributions in this paper without any parametric model assumption.
We first investigate the bias caused by semi-competing risks. Let Y be the response outcomes. In the presence of semicompeting risks, Y has four levels, with Y ¼ 1, if no adverse event occurs; Y ¼ 2, if the disease progression (DP) event occurs first; Y ¼ 3, if only the dose-limiting toxicity (DLT) event occurs; Y ¼ 4, if the DLT event occurs first followed by the DP event. Let d be the dose levels under investigation, we can define
, which can be naturally estimated from the observed data.
We use Y P and Y T to denote the binary DP and DLT outcomes.
The relationships between p jk and π lsk are:
Noticing that π lsk are the real parameters of interest and p jk are the parameters subject to semi-competing risks issue. Therefore, when the issue of the semi-competing risks arises (ϕ k > 0), it will distort the estimates of the joint probabilities π 10k and π 11k . In terms of the marginal estimates, let us define q Pk ¼ π 10k þ π 11k as the DP rate and q Tk ¼ π 01k þ π 11k as the DLT rate. Then, if the semi-competing risks issue is ignored, the DP rate and DLT rate will be represented by e q Pk ¼ p 2k þ p 4k ¼ q Pk and e q Tk ¼ p 3k þ p 4k ¼ q Tk À ϕπ 11k < q Tk so the estimate for the DLT rate will also be biased. Since the toxicity outcome (DLT) is necessary for almost any phase I and phase I/II dose-finding trial, the issue of the semicompeting risks must be appropriately addressed. The semi-competing risks problem is often nested with the late-onset issue for clinical trials because of the long follow-up time of the response outcomes. We now investigate the late-onset issue. First of all, in Figure 1, we provide an illustrative example to demonstrate the late-onset mechanism. In Figure 1, a new patient will be enrolled in the trial at the beginning of every month. Once enrolled, that patient will be followed for 2 months to ascertain the DP and DLT information. When the second patient enters the trial, the first patient has only been followed for 1 month and does not experience DP or DLT, so the late-onset issue rises. When the third patient enters the trial, the first patient has been fully followed with Y P ¼ 0 and Y T ¼ 0. However, the second patient is still under follow-up and is subject to the lateonset issue. The third patient's toxicity-efficacy information is missing when the last patient comes in due to the lateonset issue. The second patient has now been fully followed and has developed DP in his/her 1.5 months follow-up time. However, his/her DLT status is missing, but that is due to the semi-competing risks issue rather than the lateonset, and that information will always be unobservable even if the second patient has been fully followed. Because the semi-competing risks can distort the estimates and the late-onset issue is also non-ignorable, 9,21 a new method is needed to solve these entangled problems, which will be developed in the following content.
Let U T and U P be the fixed total follow-up time for DLT and DP. We define v Tki and v Pki as the actual follow-up time for DLT and DP for the ith patient assigned to dose level k at an interim analysis time and we have v Tki ≤ U T and v Pki ≤ U P by trial design. Let X Tki and X Pki be the underlying true DLT occurring time and true DP occurring time for the ith patient assigned to dose k.
, to circumvent the missing data problem caused by semicompeting risks and late-onset, we propose to approximate the probability Pr X Tki > v Tki ð Þas:
This approximation assumes that the time-to-toxicity outcome is uniformly distributed over the entire follow-up period 0, U T ½ . According to previous studies, 22 this uniform scheme is remarkably robust and yields desirable operating
, along the same line, we further approximate the probabilities Pr
After some algebra, we have the following expressions:
where IðÞ is the indicator function. We denote
Þ , 0 otherwise. Hence, the complete likelihood function in the presence of late-onset semi-competing risks outcomes can be expressed as
We develop a data augmentation method to derive the posterior distributions of the parameters q Tk , τ 0k and τ 1k based on the likelihood function (4). We assign q Tk , τ 0k and τ 1k the following Beta prior distributions:
Then, we use ω ki to denote the missing value of I X Tki > U T ð Þ . The data augmentation method is summarized as follows:
1. Impute ω ki from the Bernoulli distribution as
2. Given all the imputed data, sequentially draw the posterior samples of q Tk , τ 0k and τ 1k from the Beta distributions as
The derivation of the data augmentation method is given in Appendix. With q Tk , τ 0k and τ 1k , we can easily calculate the posterior distributions of the joint probability π lsk , which will be used in the following dose-finding design.
The goal of a phase I/II trial is to identify the OBD, which yields the best risk-benefit tradeoff among all the candidate doses. There are various ways to define the OBD (e.g., marginal probabilities, tradeoff contour, utility function), and the utility function is arguably the most convenient and popular one, which maps the multidimensional outcomes into a single index to measure the desirability of a dose in terms of risk-benefit tradeoff. The utility function has been used in a variety of phase I/II designs to define the OBD. [15][16][17]20,21,[23][24][25][26] Moreover, the utility function-based definition is flexible and incorporates other definitions as special cases. 16 We use the utility function-based definition in this paper.
To construct the utility function we need to elicit from physicians the utility values ζ ls associated with the probabilities π lsk for l, s ¼ 0, 1. The requirement for the utility values is ζ 11 < ζ 01 , ζ 10 < ζ 00 because no adverse event happening is most desirable and experiencing both DLT and DP is least desirable. Constructing the utilities requires close collaboration between clinicians and statisticians and should be customized to reflect the clinical needs and practice best. In our experience, the process of elicitation of utility values is quite natural. It can be done by simply explaining what the utilities represent to the clinicians during the decision process and asking them to specify all values after fixing the ones for the worst and best outcome. In particular, we can first fix two boundaries as ζ 11 ¼ 0 and ζ 00 ¼ 100 , and then elicit the utility values for the remaining ones by using the two boundaries as a reference, which must be located between 0 and 100. After eliciting the utility values, the utility function for dose level d ¼ k can be constructed as
ζ ls π lsk :
During the trial, given the interim data D, we can calculate the posterior mean utility as
In addition to e U kjD ð Þ, we also need to construct an admissible set to safeguard patients from overly toxic and less efficacious doses. We define ε T as the maximum tolerable DLT rate and ε P as the maximum tolerable DP rate. We first define the admissible set of toxicity as
If any dose k is found to be overly toxic during the interim analysis, then all the doses no less than k should be excluded from the trial because the toxicity rate generally increases as the dose level increases. Then, Let q Pk ¼ Pr Y P ¼ 1jd ¼ k ð Þ , we define the admissible set as
C T and C P are cut-off values, which can be calibrated through simulation studies. Finally, the dose that yields the highest posterior mean utility e U kjD ð Þ within the admissible set A is the identified OBD. Based on the aforementioned content, the toxicity rate q Tk will be underestimated, and risky doses may be included in the admissible set if the issue of the semi-competing risks is ignored. In terms of the utility function, the impact of semi-competing risks is summarized in the following theorem.
Theorem 1. The utility function U kjπ k ð Þ is overestimated if the issue of the semi-competing risks is ignored.
The proof of Theorem 1 is given in the Appendix. As the semi-competing risks can distort both A and U kjπ k ð Þ, the probability model and data augmentation algorithm developed in this paper are necessary.
We propose the following dose-finding algorithm to find the OBD for a phase I/II trial with the utility function and admissible set constructed.
1. The first cohort of patients is treated at the lowest dose level or another physician-specified dose level. 2. Update the posterior distribution of all the parameters of interest at the current dose level. 3. Use the updated posterior distribution to construct the admissible set and identify the OBD. If the admissible set is empty, early terminate the trial and conclude that no dose can be selected as the OBD. 4. Treat the next cohort of patients at the identified OBD. 5. Repeat steps 1-4 until the maximum sample size is reached. Dose skipping is not allowed when a dose-escalation occurs.
The proposed probability model does not model the untried dose. As a result, the dose-finding procedure is restricted to tired doses and may miss the global optimal dose. To solve this problem, we propose an additional dose-escalation rule. Specifically, let k à be the highest dose that has been tired so far. If during any dose-finding step, we have k à A T , indicating that the true OBD may be greater than k à , we will do the dose-escalation and allocate the next cohort of patients at dose level k à þ 1 if it is possible to further explore the dose-response profiles. Once the trial ends, we select the dose with the highest posterior mean utility within the admissible set as the recommended OBD. As an essential step to apply the proposed design to the motivating AML trial, we first evaluate the operating characteristics of the proposed design. Reporting operating characteristics is often required in trial protocols when a new design is involved. We considered a phase I/II trial with five doses and a maximum sample size of 60 in cohorts of size 3. We specified a DP upper bound ε P ¼ 0:5 with a cut-off value C P ¼ 0:9, a DLT upper bound of ε T ¼ 0:3 with a cut-off value C T ¼ 0:95. We elicited the utility values as ζ 00 ¼ 100, ζ 01 ¼ 50, ζ 10 ¼ 25 and ζ 11 ¼ 0. We specify independent noninformative beta priors Beta 0:5, 0:5 ð Þ for q Tk , τ 0k and τ 1k . Table 1 shows the simulation results with different doseresponse curves, including the OBD selection probability, the average number of patients treated at each dose level, the F I G U R E 2 Sensitivity analysis with different time-to-event generating functions average percentages of patients experiencing DP, DLT, the average sample size N and the average trial duration across 5000 simulated trials. A boldface font emphasizes the result for the true OBD. We specify a fixed total follow-up time of 2 months for every patient and an inter-arrival time τ ¼ 1 month. A new cohort of patients will enter the trial every month. We first simulate the time-to-DLT outcome from a Weibull distribution to generate the survival outcomes. Then, we simulate two independent time-to-DP outcomes conditional on Y T ¼ 0 or 1, also from the Weibull distributions. The parameters in the Weibull distributions are calibrated to match the pre-specified response rates. We compare the proposed SCI design with two alternatives. The first one uses the observable data only during the interim analysis, and we name it the observed-data design. Any patients with missing DP or DLT are ignored under the observed-data design. The second one suspends the trial during the interim analysis in the presence of late-onset outcome and reopens the trial until all the patients have been completely followed, and we name it the complete-data design. In scenario 1, none of the doses is acceptable due to overly toxic or less efficacious. All the designs under comparison yield almost surely probabilities of terminating the trial early with no dose selected. When the true OBD does exist (scenarios 2-6), as we expect, the complete-data design always yields the best performances in terms of OBD selection and patient allocation because of the fully observable data being used. However, compared with the SCI and observeddata designs, the complete-data design doubles the trial duration when the true OBD exists, making the design impractical. It is also highly unethical to frequently suspend the trial and let the patients already in the trial wait until all the data become observable. The SCI design outperforms the observed-data design across all the scenarios under consideration, and the difference can be substantial. For example, under scenario 2, the true OBD selection rate for the SCI design is 49.3%, which is about 18% higher than that for the observed-data design. Also, the number of patients at the OBD for the SCI design is 8.5 higher than the number for the observed-data design. The results for scenarios 3-6 are similar.
We also conducted two sensitivity analyses to investigate the robustness of the proposed SCI method. In Figure 2, we report the OBD selection percentages of the observed-data, SCI, and complete-data designs with different data generating distributions such as the Weibull, Lognormal, and Gamma. Figure 3 reports the OBD selection percentages of all three designs with different inter-arrival times τ. The results show that the SCI performs robustly under different data generating distributions and inter-arrival time τ. In conclusion, the pattern of the design comparison keeps consistent across all the tables and figures.
We develop a phase I/II clinical trial design to identify the OBD of IT in the presence of semi-competing risks. We propose reconstructing the likelihood based on the actual follow-up time and developing a data augmentation method to derive the posterior distributions for the parameters of interest. Simulation results confirm the desirable performance of the proposed SCI design. The complete-data design significantly prolongs the trial duration and should not be considered. Considering that the proposed SCI design is more complicated than the naive observeddata design, the latter can be considered as an alternative to the SCI design in the presence of a mild late-onset issue.
The issue of semi-competing risks has been intensively studied in the statistical community, and a lot of sophisticated survival analysis models have been proposed, such as the multi-state model 27,28 and latent failure times model, 8,[29][30][31] which treat the data as time-to-event. The benefits of treating the semi-competing risks data as time-to-event and using the survival modeling approaches include incorporating the length of time into consideration, handling censored data, and estimating the hazard ratio for association studies, among others. However, the high learning cost impedes the application of these complicated models in practical trial conduct. Moreover, the performances of the survival models rely on a series of parametric model assumptions, which are hard to justify and verify for a dose-finding trial with a limited sample size. On the other hand, the SCI design treats the censored data as missing and imputes them from the multinomial distribution for categorical data, which is more transparent and user-friendly for the clinical community. The only model assumption for the SCI design is the uniform distribution on the censored data, which is proved to be a fair and reasonable assumption. 20,22,32 Indeed, the censoring status is dynamically updated among different interim decision-making times. As the trial goes on, more and more patients change their status from censoring to observed, so the uniform distribution assumption is gradually relieved and eventually removed for the final optimal dose selection. The SCI design makes no parametric assumption on the doseresponse curve and toxicity-efficacy correlation and therefore is a robust design against parametric model assumption. Another feature of the SCI design is that it provides a unified solution for both the semi-competing and late-onset outcomes. The proposed data augmentation method treats the censoring data due to either the semi-competing risks or late-onset as the same type of missing data and imputes them using one Gibbs sampling algorithm. As a result, the SCI design can be directly applied for a mixture of semi-competing risks and/or late-onset data.
One possible extension of the SCI design is to incorporate the immune response. Immune response measures IT's biological efficacy in activating the immune system and is closely correlated with the toxicity and efficacy outcomes. By accounting for the immune response, we may further improve the performance of the SCI design. Also, IT is often used together with other drugs. Therefore, it is vital to extend the SCI design to identify the optimal dose pair for the drugdrug combination trial.
Given the Beta prior Beta α q ,β q , the posterior samples of q Tk are drawn from q Tk $ Beta
When ω ki is imputed, the conditional likelihood function for τ 0k is
Similar to Lin and Yuan, 31 we perform Taylor expansion of 1
. Consequently, we approximate L τ0 as
Based on the expression of e L τ0 , τ 0k follows the conjugate beta posterior distribution as
Along the same line, τ 1k follows the conjugate beta posterior distribution as
If the issue of the semi-competing risks is ignored, the utility function is constructed as