rpact: Confirmatory Adaptive Clinical Trial Design and Analysis

Summary

This R markdown file provides an example for planning a trial with a survival endpoint using rpact thereby illustrating the different ways of entering recruitment schemes. It also demonstrates the use of the survival simulation function.

Load rpact package

# Load rpact
library(rpact)
packageVersion("rpact") # version should be version 2.0.3 or later
## [1] '2.0.5.9009'

1 The design

A survival trial is planned to be performed with one interim stage and using an O’Brien & Fleming type \(\alpha\)-spending approach at \(\alpha = 0.025\). The interim is planned to be performed after half of the necessary events were observed. It is assumed that the median survival time is 18 months in the treatment group, and 12 months in the control. We assume that the drop-out rate is 5% after 1 year and the drop-out time is exponentially distributed.

1.1 Accrual and follow-up time given

The patients should be recruited within 12 months assuming uniform accrual. We assume an additional follow-up time of 12 months, i.e., the study should be conducted within 2 years. We also assume the survival time to be exponentially distributed.

In this simplest example, accrual and follow-up time needs to be specified and the necessary number of events and patients (total and per month) in order to reach power 90% with the assumed median survival times will be calculated.

The effect size is defined in terms of lambda1 and lambda2 (you can also specify lambda2 and hazardRatio) and the function getLambdaByMedian is used (for lambda2 the direct calculation is illustrated):

dGS <- getDesignGroupSequential(kMax = 2, typeOfDesign = "asOF", beta = 0.1)

x1 <- getSampleSizeSurvival(dGS, lambda1 = getLambdaByMedian(18), lambda2 = log(2)/12,
    dropoutRate1 = 0.05, dropoutRate2 = 0.05, dropoutTime = 12,
    accrualTime = 12, followUpTime = 12)
x1  
## Design plan parameters and output for survival data:
## 
## Design parameters:
##   Significance level                           : 0.0250 
##   Type II error rate                           : 0.1 
##   Test                                         : one-sided 
## 
## User defined parameters:
##   Lambda (1)                                   : 0.0385 
##   Lambda (2)                                   : 0.0578 
##   Follow up time                               : 12.00 
##   Drop-out rate (1)                            : 0.050 
##   Drop-out rate (2)                            : 0.050 
## 
## Default parameters:
##   Type of computation                          : Schoenfeld 
##   Theta H0                                     : 1 
##   Planned allocation ratio                     : 1 
##   Event time                                   : 12 
##   Accrual time                                 : 12.00 
##   Kappa                                        : 1 
##   Piecewise survival times                     : 0.00 
##   Drop-out time                                : 12.00 
## 
## Sample size and output:
##   Direction upper                              : FALSE 
##   pi (1)                                       : 0.370 
##   pi (2)                                       : 0.500 
##   Hazard ratio                                 : 0.667 
##   Maximum number of subjects                   : 465.7 
##   Maximum number of subjects (1)               : 232.8 
##   Maximum number of subjects (2)               : 232.8 
##   Maximum number of events                     : 256.5 
##   Accrual intensity                            : 38.8 
##   Calculate follow up time                     : FALSE 
##   Information rates [1]                        : 0.500 
##   Information rates [2]                        : 1.000 
##   Analysis times [1]                           : 13.14 
##   Analysis times [2]                           : 24.00 
##   Expected study duration under H1             : 21.26 
##   Maximal study duration                       : 24.00 
##   Number of events by stage [1]                : 128.3 
##   Number of events by stage [2]                : 256.5 
##   Expected number of events under H0           : 256.3 
##   Expected number of events under H0/H1        : 252.1 
##   Expected number of events under H1           : 224.1 
##   Number of subjects [1]                       : 467.3 
##   Number of subjects [2]                       : 467.3 
##   Expected number of subjects under H1         : 467.3 
##   Reject per stage [1]                         : 0.253 
##   Reject per stage [2]                         : 0.647 
##   Early stop                                   : 0.253 
##   Critical values (effect scale) [1]           : 0.593 
##   Critical values (effect scale) [2]           : 0.782 
##   Local one-sided significance levels [1]      : 0.001525 
##   Local one-sided significance levels [2]      : 0.024500 
## 
## Legend:
##   (i): values of treatment arm i
##   [k]: values at stage k

We conclude that a maximum number of 257 events needs to be observed to reach 90% power. Under the assumed accrual and follow-up time (and the median survival times) this number of events are expected to be observed if 466 subjects are observed at the end of the trial. This corresponds to about 38.8333333 (= accrualIntensity) patient per month to be accrued.

The interim analysis is performed after 128.3 events (which in practice needs to be rounded up or down). This number is expected to be observed at analysis time 13.14.

Note that in this case you can also enter the accrual time as accrualTime = c(0,12) and so accrualTime = 12 is just a shortcut (see below).

1.2 Accrual time and maximum number of patients given

Assume that accrual stops after 16 months with 25 patients per month, i.e., after 400 patients were recruited. We now can estimate the necessary follow-up time such that 257 events will be observed at the end of the trial.

There are two ways of defining the situation of given maximum number of subjects: You can either specify accrualTime and maxNumberOfSubjects directly; or you specify accrualTime and accrualIntensity.

In both cases, at given accrual time and number of patients, the follow-up time is calculated.

This is the code and the output for the first case:

x3 <- getSampleSizeSurvival(dGS, lambda1 = log(2)/18, lambda2 = log(2)/12,
    dropoutRate1 = 0.05, dropoutRate2 = 0.05, dropoutTime = 12,
    accrualTime = c(0, 16), accrualIntensity = 25)

    ceiling(x3$maxNumberOfSubjects)
## [1] 400
    x3$followUpTime
## [1] 15.96226
    x3$analysisTime
##          [,1]
## [1,] 16.82864
## [2,] 31.96226

This is the code and the output for the second case:

x3 <- getSampleSizeSurvival(dGS, lambda1 = log(2)/18, lambda2 = log(2)/12,
    dropoutRate1 = 0.05, dropoutRate2 = 0.05, dropoutTime = 12,
    accrualTime = c(0, 16), maxNumberOfSubjects = 400)

    x3$accrualIntensity
## [1] 25
    x3$followUpTime
## [1] 15.96226
    x3$analysisTime
##          [,1]
## [1,] 16.82864
## [2,] 31.96226

Note that for these two cases entering the accrual time in the form of c(0,16) is mandatory, i.e., you cannot enter just the end of accrual which was the shortcut before.

Alternatively, you might enter the accrualTime as a list in the form

atList <- list(
    "0 - < 16"   = 25)

and then enter accrualTime = at in getSampleSizeSurvival():

x3 <- getSampleSizeSurvival(dGS, lambda1 = log(2)/18, lambda2 = log(2)/12,
    dropoutRate1 = 0.05, dropoutRate2 = 0.05, dropoutTime = 12,
    accrualTime = atList)

    x3$accrualIntensity
## [1] 25
    x3$followUpTime
## [1] 15.96226
    x3$analysisTime
##          [,1]
## [1,] 16.82864
## [2,] 31.96226

1.3 Follow-up time and absolute intensity given

Assume now that 25 patients can be recruited each month and that there is uniform accrual. We now want to estimate the necessary accrual time if the planned follow-up time is given and equal to 12 as before.

Here the end of accrual and the number of patients is calculated at given follow-up time and absolute accrual intensity:

x2 <- getSampleSizeSurvival(dGS, hazardRatio = 2/3, lambda2 = log(2)/12,
    dropoutRate1 = 0.05, dropoutRate2 = 0.05, dropoutTime = 12,
    accrualTime = 0, accrualIntensity = 25, followUpTime = 12)

ceiling(x2$maxNumberOfSubjects)
## [1] 435
x2$accrualTime
## [1] 17.38334
x2$analysisTime
##          [,1]
## [1,] 16.79806
## [2,] 29.38334

1.4 Staggered patient entry

This case illustrates how the accrual of subjects can be entered if there is piecewise accrual over predefined intervals where it is assumed that accrual within the intervals is uniform.

For example, assume that in the first 3 months 15 patients, in the second 3 months 20 patients, and after 6 months 25 patients per month can be accrued.

If the follow-up time at given number of patients is to be calculated (Section 1.2), it is simply to enter accrualTime = c(0, 3, 6, 16) and accrualIntensity = c(15, 20, 25) as follows (maxNumberOfSubjects will be calculated):

x5 <- getSampleSizeSurvival(dGS, lambda1 = log(2)/18, lambda2 = log(2)/12,
    dropoutRate1 = 0.05, dropoutRate2 = 0.05, dropoutTime = 12,
    accrualTime = c(0, 3, 6, 16), accrualIntensity = c(15, 20, 25))

ceiling(x5$maxNumberOfSubjects)
## [1] 355
x5$followUpTime
## [1] 23.60973
x5$analysisTime
##          [,1]
## [1,] 18.83368
## [2,] 39.60973

Alternatively, you can enter

atList <- list(
    "0 - < 3" = 15,
    "3 - < 6" = 20,
    "6 - <= 16" = 25)

and run

x5 <- getSampleSizeSurvival(dGS, lambda1 = log(2)/18, lambda2 = log(2)/12,
    dropoutRate1 = 0.05, dropoutRate2 = 0.05, dropoutTime = 12,
    accrualTime = atList)

ceiling(x5$maxNumberOfSubjects)
x5$followUpTime
x5$analysisTime

If maxNumberOfSubjects is specified, the corresponding list definition is as follows:

atList <- list(
    "0 - < 3" = 15,
    "3 - < 6" = 20,
    "6 - Inf" = 25)

and

x5 <- getSampleSizeSurvival(dGS, lambda1 = log(2)/18, lambda2 = log(2)/12,
    dropoutRate1 = 0.05, dropoutRate2 = 0.05, dropoutTime = 12,
    accrualTime = atList, maxNumberOfSubjects = 355)

ceiling(x5$accrualTime)
x5$followUpTime
x5$analysisTime

does the job.

If the follow-up time is given and the end of accrual it to be calculated (Section 1.3), this can be achieved as follows:

x4 <- getSampleSizeSurvival(lambda1 = log(2)/18, lambda2 = log(2)/12,
    accrualTime = c(0, 3, 6), accrualIntensity = c(15, 20, 25), followUpTime = 12)
    
ceiling(x4$maxNumberOfSubjects)
## [1] 329
x4$accrualTime
## [1]  3.00000  6.00000 14.94986
x4$analysisTime
##          [,1]
## [1,] 26.94986

2 Verify results by simulation

Assume that the study is planned with 257 events and 400 patients under the assumptions that accrual stops after 16 months with 25 patients per month and the staggered patient entry from above. We now want to verify by simulation the correctness of the results obtained by the analytical formulae.

For this, we first use the function getPowerSurvival() in order to obtain the test characteristics if the study is performed with 257 events. The analysis time is obtained by the analytical formulae and we verify that the power is indeed exceeding 90% (note that we have to specify directionupper = FALSE in order to indicate that p-values are getting small for negative values of the test statistics, i.e., for observed hazard ratios smaller than 1):

y3 <- getPowerSurvival(dGS, lambda1 = log(2)/18, lambda2 = log(2)/12,
    dropoutRate1 = 0.05, dropoutRate2 = 0.05, dropoutTime = 12,
    accrualTime = c(0, 3, 6, 16), accrualIntensity = c(15, 20, 25),
    maxNumberOfEvents = 257, directionUpper = FALSE)

y3$analysisTime
##          [,1]
## [1,] 18.85699
## [2,] 39.74904
y3$overallReject  
## [1] 0.9005251

Practically the same result is obtained with the simulation function getSimulationSurvival() (note that plannedEvents = c(129,257) needs to be specified in teh simulation tool instead of maxNumberOfEvents = 257 in the power calculation too:

maxiter <- 10000

z3 <- getSimulationSurvival(dGS, lambda1 = log(2)/18, lambda2 = log(2)/12,
    dropoutRate1 = 0.05, dropoutRate2 = 0.05, dropoutTime = 12, maxNumberOfIterations = maxiter,
    accrualTime = c(0, 3, 6, 16), accrualIntensity = c(15, 20, 25),
    plannedEvents = c(129, 257), directionUpper = FALSE, seed = 12345)

z3$analysisTime
##          [,1]
## [1,] 18.91939
## [2,] 39.63640
z3$overallReject  
## [1] 0.9054

3 Assess sample size reassessment in adaptive survival design

Assume now that a sample size increase up to a ten-fold of the originally planned number of events is foreseen. Conditional power 90% based on the observed hazard ratios is used to increase the number of events. We can now assess by simulation the magnitute of power increase when using the appropriate method.

Note that we have to make sure that enough subjects are used in the simulation. For this, we set maxNumberOfSubjects = 3000 and no drop-outs.

Since we perform a data driven sample size recalculation (SSR), we define an inverse normal design with the same parameters as the original group sequential design:

dIN <- getDesignInverseNormal(kMax = 2, typeOfDesign = "asOF", beta = 0.1)

The following code simulates the SSR with conditional power 90% under the observed hazard ratio. minNumberOfEventsPerStage = c(NA,128), and maxNumberOfEventsPerStage = 10*c(NA,128) defines:

  • no sample size decrease

  • up to a 10 fold increase of number of events for the final stage.

This is compared with the situation without SSR:

#  Without SSR
z0 <- getPowerSurvival(dIN, lambda1 = log(2)/18, lambda2 = log(2)/12,
    maxNumberOfIterations = maxiter,
    accrualTime = c(0,16), maxNumberOfSubjects = 3000, maxNumberOfEvents = 257, 
    directionUpper = FALSE)
## Warning: Argument unknown in getPowerSurvival(...): 'maxNumberOfIterations'
## = 10000 will be ignored
z0$analysisTime
##          [,1]
## [1,] 5.584142
## [2,] 8.051325
z0$overallReject  
## [1] 0.9005251
#  With SSR
z4 <- getSimulationSurvival(dIN, lambda1 = log(2)/18, lambda2 = log(2)/12,
    maxNumberOfIterations = maxiter,
    accrualTime = c(0,16), maxNumberOfSubjects = 3000, plannedEvents = c(129, 257), 
    directionUpper = FALSE, conditionalPower = 0.9, 
    minNumberOfEventsPerStage = c(NA,128), 
    maxNumberOfEventsPerStage = 10 * c(NA,128), seed = 23456)

z4$analysisTime
##           [,1]
## [1,]  5.590061
## [2,] 11.375353
z4$overallReject  
## [1] 0.9867

Finally, we simulate the Type I error rate when using

  • the group sequential method

  • the inverse normal method

for analysing the data. The difference here is that the group sequential design uses the logrank tests for the test decision whereas with the combination test the increments of the sequence of logrank statistics are combined with the use of the inverse normal combination function. This assures Type I error control also the data-driven SSR.

The following simulation compares the Type I error rate of the inverse normal method with the Type I error rate of the (inappropriate) use of the group sequential method:

maxiter <- 10000

dGS <- getDesignGroupSequential(kMax = 2, typeOfDesign = "asOF")
dIN <- getDesignInverseNormal(kMax = 2, typeOfDesign = "asOF")
    
IN <- getSimulationSurvival(dIN, hazardRatio = 1,
    maxNumberOfIterations = maxiter,
    accrualTime = c(0, 16), maxNumberOfSubjects = 3000, plannedEvents = c(129, 257), 
    directionUpper = FALSE, conditionalPower = 0.9, 
    minNumberOfEventsPerStage = c(NA, 128), 
    maxNumberOfEventsPerStage = 10 * c(NA,128), seed = 34567)
    
GS <- getSimulationSurvival(dGS, hazardRatio = 1,
    maxNumberOfIterations = maxiter,
    accrualTime = c(0,16), maxNumberOfSubjects = 2000, plannedEvents = c(129, 257), 
    directionUpper = FALSE, conditionalPower = 0.9, minNumberOfEventsPerStage = c(NA,128), 
    maxNumberOfEventsPerStage = 10 * c(NA,128), seed = 45678)

IN$overallReject  
## [1] 0.026
GS$overallReject  
## [1] 0.0352

We see that indeed the Type I error rate is controlled for the inverse normal but nor for the group sequential method (see also this vignette).

4 References

  1. Gernot Wassmer and Werner Brannath, Group Sequential and Confirmatory Adaptive Designs in Clinical Trials, Springer 2016, ISBN 978-3319325606

System: rpact 2.0.5.9009, R version 3.6.1 (2019-07-05), platform: x86_64-w64-mingw32

To cite package ‘rpact’ in publications use:

Gernot Wassmer and Friedrich Pahlke (2019). rpact: Confirmatory Adaptive Clinical Trial Design and Analysis. R package version 2.0.5.9009. https://www.rpact.org

 

Creative Commons License
This work by Gernot Wassmer and Friedrich Pahlke is licensed under a Creative Commons Attribution-ShareAlike 4.0 International License.