OneArm2stage

Introduction

The OneArm2stage package serves multiple purposes, including:

Planning sample sizes and other design parameters for single-arm phase II trials with time-to-event (TOE) endpoints.
Conducting interim analyses, which involve making Go/No-go decisions, as well as the final analysis in survival trials.

Clinical researchers can leverage the functions within this package when conducting clinical trials focused on immunotherapies or molecularly targeted therapies for cancer, where the TOE endpoint is commonly chosen as the primary outcome measure. One notable advantage of this package is its ability to calculate design power with increased precision by utilizing the exact variance in the one-sample log-rank test (Wu et al., 2020). Furthermore, it offers multiple design solutions simultaneously for single-arm phase II survival trials, including optimal, minmax, and admissible two-stage designs (Wu et al., 2020 and Jung, 2004).

The optimal method is designed to minimize the expected sample size (ES) under H0, EN = PS_H₀ × n₁ + (1 − PS_H₀) × n, where PS_H₀ is the early stopping probability under the null hypothesis, n₁ and n are the sample size for the first stage and total sample size, respectively. In contrast, the minmax method minimizes the total sample size n.

The admissible methods aim to find a middle ground between the optimal and minmax methods by minimizing the following loss function: $$\rho=q\times n_d + (1-q) \times EN_d \space. $$ In this equation, q represents the relative importance between the maximum sample size (n) and the expected sample size under H0 (ES). The variable d signifies each individual design, with its total sample size ranging from that of the minmax method to the optimal method under the same trial conditions.

The admissible method within this package suggests selecting the design with the lowest ρ value from the options considered (refer to Example 1). This method proves especially useful in scenarios where a significant reduction in either the expected sample size (ES) or the total sample size (n) can be accomplished without significantly compromising the other parameter. For additional information regarding the admissible method, please consult Jung, 2004.

The selection among the three design methods mentioned can be determined by considering the distinct characteristics of each trial, taking into account both statistical considerations and clinical/administrative requirements. Additionally, the package allows for both restricted and unrestricted follow-up. For instance, in a design with a restricted follow-up of 3 years, each patient will be followed for 3 years (provided there are no events within the 3 years). Furthermore, any patient who reaches the end of this follow-up period without experiencing an event will be censored at the end of the 3 years for the interim and final analyses. In contrast, in a design with unrestricted follow-up also set at 3 years, the requirement is for the last patient enrolled in the trial to be followed for 3 years (provided there are no events). However, any patient enrolled before the last one will not be censored for the interim and final analyses if they have been followed for more than 3 years.

The following vignette illustrates the process of designing single-arm two-stage survival trials and conducting analyses/simulations using this package.

Examples

library(OneArm2stage)
library(survival)
library(IPDfromKM)

Example 1 (Designing a Study with Unrestricted Follow-Up)

Let’s consider the scenario where we’re designing a study with the primary endpoint being 3-year event-free survival (EFS) for sample size planning. The hypotheses under consideration are as follows:

Null Hypothesis (H0): S0 = 72%
Alternative Hypothesis (H1): S1 = 86%

These hypotheses pertain to the landmark time point of 3 years (x0=3). This is equivalent to assuming a hazard ratio of 0.459, provided that survival time follows the proportional hazard assumption (0.72^0.459 = 0.86). In this particular study, we assume a historical accrual rate of 20 patients per year, and each patient has an unrestricted follow-up time of 3 years (restricted=0 and tf=3).

Given a type I error rate of 0.05 and a desired power of 0.8, with failure times following a Weibull distribution with shape = 1, we can employ the function phase2.TTE to plan this study, as demonstrated below.

Design1 <- phase2.TTE(shape=1,S0=0.72,x0=3,hr=0.459,tf=3,rate=20,alpha=0.05,
                     beta=0.2, restricted=0)

Design1
# $Two_stage_Optimal
#   n1      c1  n      c     t1 MTSL      ES     PS
# 1 22 -0.7215 46 1.5952 1.0619  5.3 40.1736 0.2353
# 
# $Two_stage_minmax
#   n1     c1  n      c     t1 MTSL      ES     PS
# 1 34 -0.731 44 1.6257 1.6997  5.2 41.6745 0.2324
# 
# $Two_stage_Admissible
#     n1      c1  n      c     t1 MTSL      ES     PS      Rho
# 88  26 -0.7274 45 1.6061 1.2512 5.25 40.3361 0.2335 42.66805
# 98  29 -1.0069 44 1.6293 1.4269 5.20 41.5723 0.1570 42.78615
# 108 21 -0.7393 46 1.5952 1.0326 5.30 40.1736 0.2298 43.08680

The output from phase2.TTE provides design solutions based on three methods. If we opt for the optimal method, we need to plan for an accrual time of 1.06 years (equivalent to enrolling 22 patients assuming an accrual rate of 20) in the 1st stage before conducting the interim analysis. If the decision is to proceed (“Go”), we will enroll an additional 24 patients (46 − 22 = 24) for the final analysis. Critical values for the 1st and 2nd stages are -0.7215 and 1.5952, respectively. Overall, this study will take a maximum of 5.3 years to complete and is expected to have a sample size of 40.1736 under H0. The probability of early stopping under H0 is 23.53%.

Similarly, if we choose the minmax method, we need to plan for an accrual time of 1.70 years (equivalent to enrolling 34 patients assuming an accrual rate of 20) in the 1st stage before conducting the interim analysis. If the decision is to proceed (“Go”), we will enroll an additional 10 patients (44 − 34 = 10) for the final analysis. Critical values for the 1st and 2nd stages are -0.731 and 1.6257, respectively. Overall, this study will take a maximum of 5.2 years to complete and is expected to have a sample size of 41.6745 under H0. The probability of early stopping under H0 is 23.24%.

The admissible method recommends the design with the smallest ρ value among the considered designs (the first row in outputs, ρ = 42.66805). If we opt for this admissible design, we need to plan for an accrual time of 1.25 years (equivalent to enrolling 26 patients assuming an accrual rate of 20) in the 1st stage before conducting the interim analysis. If the decision is to proceed (“Go”), we will enroll an additional 19 patients (45 − 26 = 19) for the final analysis. Critical values for the 1st and 2nd stages are -0.7274 and 1.6061, respectively. Overall, this study will take a maximum of 5.25 years to complete and is expected to have a sample size of 40.3361 under H0. The probability of early stopping under H0 is 23.35%.

Example 2 (When Early Stopping Probability is Pre-Defined)

If there is a need to pre-define the early stopping probability under H0, you can do so by using the argument prStop. This argument defines the the lower limit of the early stopping probability under H0. The default value for “prStop” is 0, which technically means the early stopping probability is > = 0. Essentially, it also means there is no restriction on the early stopping probability as any probability should be equal or larger than 0. In the following example, we set the early stopping probability to be at least 35%.

Design2 <- phase2.TTE(shape=1,S0=0.72,x0=3,hr=0.459,tf=3,rate=20,alpha=0.05,
                     beta=0.2, prStop=0.35, restricted=0)

Design2
# $Two_stage_Optimal
#   n1      c1  n      c     t1 MTSL      ES     PS
# 1 29 -0.3844 47 1.5799 1.4364 5.35 40.5985 0.3503
# 
# $Two_stage_minmax
#   n1      c1  n      c     t1 MTSL      ES     PS
# 1 40 -0.3829 44 1.6171 1.9808  5.2 42.4619 0.3509
# 
# $Two_stage_Admissible
#     n1      c1  n      c     t1 MTSL      ES     PS      Rho
# 92  34 -0.3850 45 1.5988 1.6732 5.25 40.9611 0.3501 42.98055
# 102 40 -0.3846 44 1.6171 1.9791 5.20 42.4522 0.3503 43.22610
# 921 31 -0.3851 46 1.5884 1.5339 5.30 40.6365 0.3501 43.31825
# 91  29 -0.3850 47 1.5799 1.4353 5.35 40.5949 0.3501 43.79745

We can compare the outputs above with those in Example 1. In each of the designs, the early stopping probability (PS) under H0 is now set at 35%. This strategy allows for a higher probability of early termination of an ineffective drug during the trial process. We should also be aware that using this approach necessitates larger sample sizes, particularly for the 1st stage. For example, in the optimal design, setting PS to 35% results in sample sizes of 29 for the 1st stage and 47 for the 2nd stage, whereas in the first example, the sample sizes were 22 for the 1st stage and 46 for the 2nd stage.

Example 3 (Exploring the Admissible Method Further)

The design derived from the admissible method can be fine-tuned by adjusting q_value, the parameter determines the relative importance assigned to the maximum sample size (n) versus the expected sample size under H0 (ES) when creating the admissible method-based design. Its default value is 0.5, with a range of (0, 1). Smaller values of q_value emphasize ES, while larger values prioritize the maximum sample size, n. Let’s compare the results using the code below.

Design3 <- phase2.TTE(shape=1,S0=0.72,x0=3,hr=0.459,tf=3,rate=20,alpha=0.05,
                     beta=0.2, q_value=0.1, restricted=0)

Design3_2 <- phase2.TTE(shape=1,S0=0.72,x0=3,hr=0.459,tf=3,rate=20,alpha=0.05,
                     beta=0.2, q_value=0.7, restricted=0)


# q_value = 0.1
Design3
# 
# $Two_stage_Admissible
#     n1      c1  n      c     t1 MTSL      ES     PS      Rho
# 108 21 -0.7393 46 1.5952 1.0326 5.30 40.1736 0.2298 40.75624
# 88  26 -0.7274 45 1.6061 1.2512 5.25 40.3361 0.2335 40.80249
# 98  29 -1.0069 44 1.6293 1.4269 5.20 41.5723 0.1570 41.81507

# q_value = 0.7
Design3_2
# $Two_stage_Admissible
#     n1      c1  n      c     t1 MTSL      ES     PS      Rho
# 98  29 -1.0069 44 1.6293 1.4269 5.20 41.5723 0.1570 43.27169
# 88  26 -0.7274 45 1.6061 1.2512 5.25 40.3361 0.2335 43.60083
# 108 21 -0.7393 46 1.5952 1.0326 5.30 40.1736 0.2298 44.25208

Modifications to the q_value solely impact the admissible method. Consequently, the designs obtained through the optimal and minmax methods remain unchanged from those presented in Example 1 (details not shown here). Compared to Design1 in Example 1: when q_value equals 0.1, the admissible design (n=46) closely resembles the design stemming from the optimal method, which prioritizes ES. Conversely, when q_value is set to 0.7, the design (n=44) more closely aligns with the design derived from the minmax method, emphasizing n.

Example 4 (Conducting Interim/Final Analysis)

Let’s assume we opt for the design based on the admissible method. We would need to plan for an accrual time of 1.2512 years (equivalent to enrolling 26 patients, assuming a rate of 20 patients per year) in the 1st stage before conducting the interim analysis. If the decision is to proceed (“Go”), we will enroll an additional 19 patients in the 2nd stage (45 – 26 = 19) to perform the final analysis. The critical values for the 1st and 2nd stages are -0.7274 and 1.6061, respectively. In total, this study will take a maximum of 5.25 years to complete, with an expected sample size of 40.3361 under H0. The probability of early termination under H0 is 23.35%.

Suppose we have adhered to the aforementioned design and gathered data from 26 patients for four variables:

Entry (the point at which patients enter the trial),
Time (the observed time that patients spend in the trial),
Status (an indicator of patient status), and
Total (the sum of Entry and Time), as illustrated below.

For those interested in trying the example themselves, the data samples used can be located in the folder at /OneArm2stage/extdata.

#>         Entry      time status     Total
#> 16 0.00140574 1.2497943      0 5.2500000
#> 5  0.10312013 1.1480799      0 5.2500000
#> 11 0.21339149 1.0378085      0 5.2500000
#> 24 0.23144545 1.0197546      0 5.2500000
#> 21 0.25005470 1.0011453      0 5.2500000
#> 7  0.27427333 0.9769267      0 5.2500000
#> 9  0.28694621 0.9642538      0 5.2500000
#> 1  0.31231364 0.9388864      0 5.2500000
#> 26 0.39386846 0.8573315      0 5.2500000
#> 8  0.46469563 0.4029726      1 0.8676682
#> 18 0.49526749 0.3440752      1 0.8393427
#> 2  0.52432672 0.7268733      0 2.6966289
#> 22 0.54814381 0.7030562      0 5.2500000
#> 4  0.59843844 0.6527616      0 1.7345921
#> 13 0.61736320 0.6338368      0 5.2500000
#> 20 0.79154530 0.4596547      0 2.8812490
#> 10 0.84254125 0.4086588      0 5.2500000
#> 19 0.85458721 0.3966128      0 1.7934731
#> 12 0.86393169 0.3872683      0 5.2500000
#> 23 0.93970525 0.3114947      0 5.2500000
#> 25 0.97850867 0.2726913      0 5.2500000
#> 15 0.98962130 0.2615787      0 4.8438482
#> 6  0.99495017 0.2562498      0 5.2500000
#> 14 1.00916177 0.2420382      0 2.0781127
#> 3  1.04841551 0.2027845      0 5.2500000
#> 17 1.06946229 0.1817377      0 5.2500000

We can create a Kaplan-Meier survival curve using this hypothetical dataset (dat1, n=26).

The survival rates, as estimated by the Kaplan-Meier method, can also be obtained using this hypothetical dataset (dat1, n=26).

#> Call: survfit(formula = Surv(time, status) ~ 1, data = dat1, conf.type = "log-log")
#> 
#>   time n.risk n.event survival std.err lower 95% CI upper 95% CI
#>  0.182     26       0    1.000  0.0000           NA           NA
#>  0.203     25       0    1.000  0.0000           NA           NA
#>  0.242     24       0    1.000  0.0000           NA           NA
#>  0.256     23       0    1.000  0.0000           NA           NA
#>  0.262     22       0    1.000  0.0000           NA           NA
#>  0.273     21       0    1.000  0.0000           NA           NA
#>  0.311     20       0    1.000  0.0000           NA           NA
#>  0.344     19       1    0.947  0.0512        0.681        0.992
#>  0.387     18       0    0.947  0.0512        0.681        0.992
#>  0.397     17       0    0.947  0.0512        0.681        0.992
#>  0.403     16       1    0.888  0.0748        0.621        0.971
#>  0.409     15       0    0.888  0.0748        0.621        0.971
#>  0.460     14       0    0.888  0.0748        0.621        0.971
#>  0.634     13       0    0.888  0.0748        0.621        0.971
#>  0.653     12       0    0.888  0.0748        0.621        0.971
#>  0.703     11       0    0.888  0.0748        0.621        0.971
#>  0.727     10       0    0.888  0.0748        0.621        0.971
#>  0.857      9       0    0.888  0.0748        0.621        0.971
#>  0.939      8       0    0.888  0.0748        0.621        0.971
#>  0.964      7       0    0.888  0.0748        0.621        0.971
#>  0.977      6       0    0.888  0.0748        0.621        0.971
#>  1.001      5       0    0.888  0.0748        0.621        0.971
#>  1.020      4       0    0.888  0.0748        0.621        0.971
#>  1.038      3       0    0.888  0.0748        0.621        0.971
#>  1.148      2       0    0.888  0.0748        0.621        0.971
#>  1.250      1       0    0.888  0.0748        0.621        0.971

We can employ the LRT (log-rank test) function in the OneArm2stage package to evaluate H0: S0 = 0.72 versus H1: S1 = 0.86 at the landmark time x0=3, utilizing the data provided earlier (n = 26). The determination of “Go” or “No-go” can be established based on the computed critical value Z from the following output.

LRT(shape = 1, S0=0.72, x0=3, data=dat1)
#>       O       E       Z 
#>  2.0000  1.7562 -0.1840

According to the output provided above, the test statistic Z is -0.1840, which is larger than the critical value for the 1st stage, c1 = -0.7274 (refer to Example 1). Hence, the decision is to proceed with a “Go,” and the trial will advance to the 2nd stage. As a result, we will enroll an additional 19 patients (45 - 26 = 19). We assume the final dataset, after adding these 19 patients, is as presented below.

#>         Entry       time status     Total
#> 29 0.00140574 5.24859426      0 5.2500000
#> 6  0.10312013 5.14687987      0 5.2500000
#> 17 0.21339149 5.03660851      0 5.2500000
#> 40 0.23144545 5.01855455      0 5.2500000
#> 35 0.25005470 4.99994530      0 5.2500000
#> 9  0.27427333 4.97572667      0 5.2500000
#> 12 0.28694621 4.96305379      0 5.2500000
#> 1  0.31231364 4.93768636      0 5.2500000
#> 45 0.39386846 4.85613154      0 5.2500000
#> 11 0.46469563 0.40297260      1 0.8676682
#> 31 0.49526749 0.34407523      1 0.8393427
#> 2  0.52432672 2.17230221      1 2.6966289
#> 36 0.54814381 4.70185619      0 5.2500000
#> 4  0.59843844 1.13615368      1 1.7345921
#> 19 0.61736320 4.63263680      0 5.2500000
#> 34 0.79154530 2.08970372      1 2.8812490
#> 15 0.84254125 4.40745875      0 5.2500000
#> 32 0.85458721 0.93888589      1 1.7934731
#> 18 0.86393169 4.38606831      0 5.2500000
#> 38 0.93970525 4.31029475      0 5.2500000
#> 41 0.97850867 4.27149133      0 5.2500000
#> 25 0.98962130 3.85422687      1 4.8438482
#> 7  0.99495017 4.25504983      0 5.2500000
#> 21 1.00916177 1.06895090      1 2.0781127
#> 3  1.04841551 4.20158449      0 5.2500000
#> 30 1.06946229 4.18053771      0 5.2500000
#> 10 1.26213296 3.98786704      0 5.2500000
#> 33 1.37873476 3.37933771      1 4.7580725
#> 27 1.41574755 3.83425245      0 5.2500000
#> 16 1.49650919 0.96396078      1 2.4604700
#> 37 1.50312507 2.52477082      1 4.0278959
#> 28 1.59791040 3.65208960      0 5.2500000
#> 13 1.69494269 3.55467264      1 5.2496153
#> 26 1.69756911 3.14664796      1 4.8442171
#> 39 1.77344063 3.47655937      0 5.2500000
#> 24 1.78727022 0.05248067      1 1.8397509
#> 8  1.79758090 1.03943440      1 2.8370153
#> 22 1.82264479 3.35163892      1 5.1742837
#> 23 1.82787640 3.42212360      0 5.2500000
#> 20 1.83294009 0.42523503      1 2.2581751
#> 5  1.93011236 0.56052654      1 2.4906389
#> 44 1.99455539 3.25544461      0 5.2500000
#> 43 2.00936501 3.24063499      0 5.2500000
#> 14 2.01385206 3.23614794      0 5.2500000
#> 42 2.21615320 3.03384680      0 5.2500000

We can construct a Kaplan-Meier survival curve for this combined dataset (n=45).

The survival rates, as estimated using the Kaplan-Meier method, can also be derived from this hypothetical dataset (dat2, n=45). It is evident that at approximately 3.03 years, the survival rate is 0.711, with a 95% confidence interval of (0.555, 0.821).

#> Call: survfit(formula = Surv(time, status) ~ 1, data = dat2, conf.type = "log-log")
#> 
#>    time n.risk n.event survival std.err lower 95% CI upper 95% CI
#>  0.0525     45       1    0.978  0.0220        0.853        0.997
#>  0.3441     44       1    0.956  0.0307        0.834        0.989
#>  0.4030     43       1    0.933  0.0372        0.807        0.978
#>  0.4252     42       1    0.911  0.0424        0.780        0.966
#>  0.5605     41       1    0.889  0.0468        0.753        0.952
#>  0.9389     40       1    0.867  0.0507        0.727        0.938
#>  0.9640     39       1    0.844  0.0540        0.701        0.923
#>  1.0394     38       1    0.822  0.0570        0.676        0.907
#>  1.0690     37       1    0.800  0.0596        0.651        0.891
#>  1.1362     36       1    0.778  0.0620        0.626        0.874
#>  2.0897     35       1    0.756  0.0641        0.602        0.856
#>  2.1723     34       1    0.733  0.0659        0.578        0.839
#>  2.5248     33       1    0.711  0.0676        0.555        0.821
#>  3.0338     32       0    0.711  0.0676        0.555        0.821
#>  3.1466     31       1    0.688  0.0692        0.531        0.802
#>  3.2361     30       0    0.688  0.0692        0.531        0.802
#>  3.2406     29       0    0.688  0.0692        0.531        0.802
#>  3.2554     28       0    0.688  0.0692        0.531        0.802
#>  3.3516     27       1    0.663  0.0712        0.503        0.781
#>  3.3793     26       1    0.637  0.0728        0.477        0.760
#>  3.4221     25       0    0.637  0.0728        0.477        0.760
#>  3.4766     24       0    0.637  0.0728        0.477        0.760
#>  3.5547     23       1    0.609  0.0748        0.447        0.737
#>  3.6521     22       0    0.609  0.0748        0.447        0.737
#>  3.8343     21       0    0.609  0.0748        0.447        0.737
#>  3.8542     20       1    0.579  0.0770        0.415        0.712
#>  3.9879     19       0    0.579  0.0770        0.415        0.712
#>  4.1805     18       0    0.579  0.0770        0.415        0.712
#>  4.2016     17       0    0.579  0.0770        0.415        0.712
#>  4.2550     16       0    0.579  0.0770        0.415        0.712
#>  4.2715     15       0    0.579  0.0770        0.415        0.712
#>  4.3103     14       0    0.579  0.0770        0.415        0.712
#>  4.3861     13       0    0.579  0.0770        0.415        0.712
#>  4.4075     12       0    0.579  0.0770        0.415        0.712
#>  4.6326     11       0    0.579  0.0770        0.415        0.712
#>  4.7019     10       0    0.579  0.0770        0.415        0.712
#>  4.8561      9       0    0.579  0.0770        0.415        0.712
#>  4.9377      8       0    0.579  0.0770        0.415        0.712
#>  4.9631      7       0    0.579  0.0770        0.415        0.712
#>  4.9757      6       0    0.579  0.0770        0.415        0.712
#>  4.9999      5       0    0.579  0.0770        0.415        0.712
#>  5.0186      4       0    0.579  0.0770        0.415        0.712
#>  5.0366      3       0    0.579  0.0770        0.415        0.712
#>  5.1469      2       0    0.579  0.0770        0.415        0.712
#>  5.2486      1       0    0.579  0.0770        0.415        0.712

Let’s employ the LRT (log-rank test) to perform the final analysis using data from a total of 45 patients.

LRT(shape=1, S0=0.72, x0=3, data=dat2)
#>       O       E       Z 
#> 18.0000 16.0611 -0.4838

Based on the provided output, we observe that the test statistic for the final stage, denoted as Z, equals -0.4838. This value is smaller than the critical value for the final stage, c = 1.6061 (refer to Example 1). Consequently, we are unable to reject the null hypothesis, and as a result, we must declare the trial as futile. Additionally, the LRT output also provides information regarding the observed and expected number of events for the final stage (O = 18 and E = 16.0611).

Example 5 (Designing a Study with Restricted Follow-Up)

Let’s consider designing a study similar to the one presented in Example 1, but this time with restricted follow-up limited to 3 years (tf = 3). We can continue to utilize the phase2.TTE function with the argument restricted=1.

Design5 <- phase2.TTE(shape=1,S0=0.72,x0=3,hr=0.459,tf=3,rate=20,alpha=0.05,
                     beta=0.2, restricted=1)

Design5
# $Two_stage_Optimal
#   n1      c1  n   c     t1 MTSL      ES    PS
# 1 36 -0.3665 60 1.6 1.7523    6 51.0914 0.357
# 
# $Two_stage_minmax
#   n1      c1  n      c     t1 MTSL      ES     PS
# 1 37 -0.5517 58 1.6193 1.8222  5.9 51.7359 0.2906
# 
# $Two_stage_Admissible
#     n1      c1  n      c     t1 MTSL      ES     PS      Rho
# 120 37 -0.5625 58 1.6196 1.8066 5.90 51.7260 0.2869 54.86300
# 88  36 -0.4531 59 1.6086 1.7515 5.95 51.2044 0.3252 55.10220
# 232 36 -0.3643 60 1.6000 1.7556 6.00 51.0949 0.3578 55.54745

If we opt to utilize the admissible design based on the aforementioned results, our accrual time for the 1st stage would be 1.8 years (equivalent to enrolling 37 patients at an accrual rate of 20). Should the interim analysis result in a “Go” decision, we will then enroll an additional 21 patients (58 − 37 = 21) during the 2nd stage before proceeding with the final analysis. The critical values for the 1st and 2nd stages are -0.5625 and 1.6196, respectively. In summary, this study is expected to take a maximum of 5.9 years to complete, with an anticipated sample size of 51.726 under H0. The probability of early termination under H0 is 28.69%. When compared to the unrestricted design (as in Example 1), it appears that the restricted design may require larger sample sizes and a longer completion time.

Example 6 (Simulation)

The OneArm2stage package also provides the Sim function, which allows us to calculate empirical power and type-I error through simulation based on the specified design parameters (e.g., t1, n1, n, c1, c) for any two-stage design obtained from phase2.TTE. The restricted argument can be utilized to specify whether follow-up is restricted (1) or unrestricted (0) in Sim.

As part of our examples, we will conduct simulations to determine the empirical power of the design corresponding to the results obtained from phase2.TTE in Example 1 using Sim.

Design1$param
# shape   S0    hr alpha beta rate x0 tf q_value prStop restricted
#     1 0.72 0.459  0.05  0.2   20  3  3     0.5      0          0

Design1$Two_stage_Optimal
#   n1      c1  n      c     t1 MTSL      ES     PS
# 1 22 -0.7215 46 1.5952 1.0619  5.3 40.1736 0.2353

Design1$Two_stage_minmax
#   n1     c1  n      c     t1 MTSL      ES     PS
# 1 34 -0.731 44 1.6257 1.6997  5.2 41.6745 0.2324

Design1$Two_stage_Admissible
#     n1      c1  n      c     t1 MTSL      ES     PS      Rho
# 88  26 -0.7274 45 1.6061 1.2512 5.25 40.3361 0.2335 42.66805
# 98  29 -1.0069 44 1.6293 1.4269 5.20 41.5723 0.1570 42.78615
# 108 21 -0.7393 46 1.5952 1.0326 5.30 40.1736 0.2298 43.08680

For unrestricted designs, you can omit the restricted=0 argument as in the following code, as this is the default setting.

# calculate empirical power based on optimal method, power=0.778
Sim(shape=1, S0=0.72, S1=0.72^(0.459), x0=3, tf=3, rate=20, t1=1.0619,
     c1=-0.7215, c=1.5952, n1=22, n=46, N=10000, seed=5868)
#> [1] 0.778

# calculate empirical power based on minmax method, power=0.803
Sim(shape=1, S0=0.72, S1=0.72^(0.459), x0=3, tf=3, rate=20, t1=1.6997,
     c1=-0.731, c=1.6257, n1=34, n=44, N=10000, seed=5868)
#> [1] 0.803

# calculate empirical power based on admissible method, power=0.789
Sim(shape=1, S0=0.72, S1=0.72^(0.459), x0=3, tf=3, rate=20, t1=1.2512,
     c1=-0.7274, c=1.6061, n1=26, n=45, N=10000, seed=5868)
#> [1] 0.789

However, if you intend to use Sim for simulating a restricted design, it is essential to explicitly specify restricted=1. Let’s try simulate the restricted design we obtained in Example 5.

Design5
# $Two_stage_Optimal
#   n1      c1  n   c     t1 MTSL      ES    PS
# 1 36 -0.3665 60 1.6 1.7523    6 51.0914 0.357
# 
# $Two_stage_minmax
#   n1      c1  n      c     t1 MTSL      ES     PS
# 1 37 -0.5517 58 1.6193 1.8222  5.9 51.7359 0.2906
# 
# $Two_stage_Admissible
#     n1      c1  n      c     t1 MTSL      ES     PS      Rho
# 120 37 -0.5625 58 1.6196 1.8066 5.90 51.7260 0.2869 54.86300
# 88  36 -0.4531 59 1.6086 1.7515 5.95 51.2044 0.3252 55.10220
# 232 36 -0.3643 60 1.6000 1.7556 6.00 51.0949 0.3578 55.54745

# calculate empirical power based on optimal method, power=0.799
Sim(shape=1, S0=0.72, S1=0.72^(0.459), x0=3, tf=3, rate=20, t1=1.7523,
     c1=-0.3665, c=1.6, n1=36, n=60, N=10000, restricted=1, seed=5868)
#> [1] 0.799

# calculate empirical power based on minmax method, power=0.796
Sim(shape=1, S0=0.72, S1=0.72^(0.459), x0=3, tf=3, rate=20, t1=1.8222,
     c1=-0.5517, c=1.6193, n1=37, n=58, N=10000, restricted=1, seed=5868)
#> [1] 0.796

# calculate empirical power based on admissible method, power=0.796
Sim(shape=1, S0=0.72, S1=0.72^(0.459), x0=3, tf=3, rate=20, t1=1.8066,
     c1=-0.5625, c=1.6196, n1=37, n=58, N=10000, restricted=1, seed=5868)
#> [1] 0.796

Example 7 (Fitting Historical Data)

Suppose we wish to plan a trial assuming that the failure time follows the Weibull distribution, but we lack knowledge about the parameters of the underlying distribution. In such a situation, we can utilize historical data from previous literature to make an informed estimation. The IPDfromKM package enables the extraction of individual patient data (IPD) from digitized Kaplan-Meier plots. The following example code demonstrates this process using a study conducted by Wang et al. (2018) on a group of mantle cell lymphoma patients. We have selected the KM curve depicting overall survival (Figure 2D in Wang et al., 2018) as the plot from which to extract the IPD.

The steps outlined below provide guidance on how to extract IPD data on your own, and you are encouraged to try it out.

## First, we take a screenshot of the plot and save it as a .png file. In this example it is named as "wang2018.png".

## Next, create a path directing to the .png file.
## e.g. path <-"./vignettes/wang2018.png"

## Last, type the following code in console, and follow the prompts to extract data points by clicking your mouse on the plot.
## df <- getpoints(path, 0, 24, 0, 100)  # 0,24,0,100 are the ranges of x and y-axes in the image

However, for the purposes of this vignette, we have provided a sample df dataset that can be used directly.

# a sample dataset that we already extracted from Wang et al, 2018.
df<- read.csv(system.file("extdata", "df.csv", package = "OneArm2stage"))

# risk time points
trisk <- c(0,2,4,6,8,10,12,14,16,18,20,22,24) 

# number of patients at risk at each risk time point
nrisk.radio <- c(124,120,115,110,107,104,103,95,46,18,11,8,0) 

# Preprocess the raw coordinates into an proper format for reco""nstruct IPD
pre_radio <- preprocess(dat=df, trisk=trisk, nrisk=nrisk.radio,totalpts=NULL,maxy=100)

#Reconstruct IPD
est_radio <- getIPD(prep=pre_radio,armID=0,tot.events=NULL)

We can conduct a survival analysis on the reconstructed IPD and compare the results with those reported in the original paper. The estimated survival probability at time = 12 is 0.86, with a 95% CI of (0.8, 0.92). This closely aligns with the result reported in Wang et al., 2018, which is 0.87 with a 95% CI of (0.79, 0.92).

plot_ex5 <- survreport(ipd1=est_radio$IPD,ipd2=NULL,arms=1,interval=2,s=c(0.75,0.5,0.25),showplots=F)
plot_ex5$arm1$survprob
#>             Surv     SE 0.95LCI 0.95UCI
#> time = 2  0.9677 0.0159  0.9371  0.9993
#> time = 4  0.9433 0.0208  0.9033  0.9850
#> time = 6  0.9184 0.0247  0.8712  0.9682
#> time = 8  0.8934 0.0280  0.8402  0.9499
#> time = 10 0.8681 0.0307  0.8099  0.9305
#> time = 12 0.8597 0.0316  0.8000  0.9239
#> time = 14 0.8254 0.0347  0.7602  0.8962
#> time = 16 0.7820 0.0379  0.7112  0.8599
#> time = 18 0.7650 0.0407  0.6893  0.8491
#> time = 20 0.6400 0.0914  0.4837  0.8467
#> time = 22 0.6400 0.0914  0.4837  0.8467

Now, using the reconstructed IPD, we can fit a parametric regression model assuming a Weibull distribution (the most widely used option). Before proceeding, we need to reformat the IPD data to make it compatible with the FitDat function in our package.

ipd <- est_radio$IPD
dat3 <- as.data.frame(cbind(rep(0, nrow(ipd)),ipd$time, ipd$status))
colnames(dat3) <- c("Entry", "Time", "Cens")

We apply the FitDat function to fit the historical data with the Weibull distribution. This function provides the AIC value for the fitted model and estimates the parameters. We can then use these parameter estimates to design the trial, as demonstrated in Examples 1 and 5.

modelWB<- FitDat(dat3)

#AIC
modelWB$AIC
#>  Weibull 
#> 301.7776

modelWB$parameter.estimates
#> $Weibull
#>     shape     scale 
#> 0.1133671 3.9939753

Reference

Wu, J, Chen L, Wei J, Weiss H, Chauhan A. (2020). Two-stage phase II survival trial design. Pharmaceutical Statistics. 19:214-229. https://doi.org/10.1002/pst.1983

Wang, M., Rule, S., Zinzani, P. L., Goy, A., Casasnovas, O., Smith, S. D.,…, Robak, T. (2018). Acalabrutinib in relapsed or refractory mantle cell lymphoma (ACE-LY-004): a single-arm, multicentre, phase 2 trial. The Lancet. 391(10121), 659–667. https://doi.org/10.1016/s0140-6736(17)33108-2

Jung, S.-H., (2004). Admissible two-stage designs for phase II cancer clinical trials. Statistics in Medicine. 23:561-569. https://doi.org/10.1002/sim.1600