Planning a Survival Trial

Daniel Sabanés Bové, Friedrich Pahlke

RPACT

June 9, 2026

Overview

  • This example is adapted from the corresponding vignette written by Gernot Wassmer and Friedrich Pahlke
  • A survival trial is planned to be performed with one interim stage and using an O’Brien & Fleming type \(\alpha\)-spending approach
  • We have further assumptions about the interim analysis, the median survival times in the two trial groups, as well as the dropout rate
  • We will illustrate the different ways of entering recruitment schemes, as well as simulations

Types of calculations

Starting point: the design object

  • One interim analysis is planned.
  • We use an O’Brien & Fleming type alpha-spending design.
  • One-sided type I error is \(\alpha = 0.025\).
  • Target power is \(90\%\).
library(rpact)
designGS <- getDesignGroupSequential(
    kMax = 2,
    typeOfDesign = "asOF",
    beta = 0.1
)
designGS

Design parameters and output of group sequential design

User defined parameters

  • Type of design: O’Brien & Fleming type alpha spending
  • Maximum number of stages: 2
  • Stages: 1, 2
  • Type II error rate: 0.1000

Derived from user defined parameters

  • Information rates: 0.500, 1.000
  • Futility bounds (non-binding): -Inf

Default parameters

  • Significance level: 0.0250
  • Binding futility: FALSE
  • Test: one-sided
  • Tolerance: 1e-08
  • Type of beta spending: none

Output

  • Cumulative alpha spending: 0.001525, 0.025000
  • Critical values: 2.963, 1.969
  • Stage levels (one-sided): 0.001525, 0.024500

Further assumptions are needed

  • Median survival time in treatment group: 18 months.
  • Median survival time in control group: 12 months.
  • Dropout rate: 5% after 1 year.
  • Dropout time is exponentially distributed.

Case 1: Accrual and follow-up time given

Assume:

  • Uniform accrual over 12 months.
  • Additional follow-up time of 12 months.
  • Exponential survival times.

We calculate the required number of events and subjects.

  • The effect can be entered via lambda1 and lambda2, with getLambdaByMedian() as a convenient helper.

Case 1: Accrual and follow-up time given (cont’d)

sampleSize <- designGS |>
    getSampleSizeSurvival(
        lambda1 = getLambdaByMedian(18),
        lambda2 = log(2) / 12,
        dropoutRate1 = 0.05,
        dropoutRate2 = 0.05,
        dropoutTime = 12,
        accrualTime = 12,
        followUpTime = 12
    )
sampleSize |> summary()

Sample size calculation for a survival endpoint

Sequential analysis with a maximum of 2 looks (group sequential design), one-sided overall significance level 2.5%, power 90%. The results were calculated for a two-sample logrank test, H0: hazard ratio = 1, H1: treatment lambda(1) = 0.039, control lambda(2) = 0.058, accrual time = 12, accrual intensity = 38.9, follow-up time = 12, dropout rate(1) = 0.05, dropout rate(2) = 0.05, dropout time = 12.

Stage 1 2
Planned information rate 50% 100%
Cumulative alpha spent 0.0015 0.0250
Stage levels (one-sided) 0.0015 0.0245
Efficacy boundary (z-value scale) 2.963 1.969
Efficacy boundary (t) 0.593 0.782
Cumulative power 0.2525 0.9000
Number of subjects 467.3 467.3
Expected number of subjects under H1 467.3
Cumulative number of events 128.3 256.5
Expected number of events under H1 224.1
Analysis time 13.14 24.00
Expected study duration under H1 21.26
Exit probability for efficacy (under H0) 0.0015
Exit probability for efficacy (under H1) 0.2525

Legend:

  • (t): treatment effect scale

Case 1: Accrual and follow-up time given (cont’d)

ceiling(sampleSize$maxNumberOfEvents) # how many events?
[1] 257
ceiling(sampleSize$maxNumberOfSubjects) # how many subjects?
[1] 468
ceiling(sampleSize$maxNumberOfSubjects) / 12 # accrual per month?
[1] 39
round(sampleSize$cumulativeEventsPerStage[1, 1], 1) # IA at how many events?
[1] 128.3
round(sampleSize$analysisTime[1, 1], 2) # When is the IA?
[1] 13.14

Shortcut for accrual time input

For this case, both inputs are equivalent:

  • accrualTime = 12
  • accrualTime = c(0, 12)

The scalar form is just a shortcut in this specific setup.

Digits settings

With rpact we have lots of flexibility in how the results are displayed, e.g. regarding the number of digits, rounding function, etc. This can be set globally for all outputs or for specific statistics. For example, if we want to display the number of subjects without decimals and always round up, we can set this globally for all outputs:

setOutputFormat("maxNumberOfSubjects",
    digits = 0,
    nsmall = 0, roundFunction = "ceiling"
)

Case 2: Accrual time and maximum subjects given

Now assume:

  • Accrual stops after 16 months.
  • 25 patients per month.

So around 400 subjects are accrued, and we estimate follow-up time.

Case 2a: specify accrual time and intensity

designGS |>
    getSampleSizeSurvival(
        median1 = 18,
        median2 = 12,
        dropoutRate1 = 0.05,
        dropoutRate2 = 0.05,
        dropoutTime = 12,
        accrualTime = c(0, 16), # not just 16 here!
        accrualIntensity = 25
    ) |>
    fetch(
        maxNumberOfSubjects,
        followUpTime,
        analysisTime
    )
$maxNumberOfSubjects
[1] 400

$followUpTime
[1] 15.96226

$analysisTime
         [,1]
[1,] 16.82864
[2,] 31.96226

Case 2b: specify accrual time and maximum subjects

designGS |>
    getSampleSizeSurvival(
        median1 = 18,
        median2 = 12,
        dropoutRate1 = 0.05,
        dropoutRate2 = 0.05,
        dropoutTime = 12,
        accrualTime = c(0, 16), # not just 16 here!
        maxNumberOfSubjects = 400
    ) |>
    fetch(
        accrualIntensity,
        followUpTime,
        analysisTime
    )
$accrualIntensity
[1] 25

$followUpTime
[1] 15.96226

$analysisTime
         [,1]
[1,] 16.82864
[2,] 31.96226

Alternative list notation for accrual scheme

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

designGS |>
    getSampleSizeSurvival(
        median1 = 18,
        median2 = 12,
        dropoutRate1 = 0.05,
        dropoutRate2 = 0.05,
        dropoutTime = 12,
        accrualTime = accrualTimetList
    ) |>
    fetch(
        accrualIntensity,
        followUpTime,
        analysisTime
    )
$accrualIntensity
[1] 25

$followUpTime
[1] 15.96226

$analysisTime
         [,1]
[1,] 16.82864
[2,] 31.96226

Case 3: Follow-up time and absolute intensity given

Now assume:

  • Absolute accrual intensity is 25 patients per month.
  • Follow-up time is fixed at 12 months.

We estimate end of accrual and total number of patients.

Case 3: Follow-up time and absolute intensity given (cont’d)

designGS |>
    getSampleSizeSurvival(
        hazardRatio = 2 / 3,
        median2 = 12,
        dropoutRate1 = 0.05,
        dropoutRate2 = 0.05,
        dropoutTime = 12,
        accrualTime = 0,
        accrualIntensity = 25,
        followUpTime = 12
    ) |>
    fetch(
        maxNumberOfSubjects,
        followUpTime,
        analysisTime
    )
$maxNumberOfSubjects
[1] 434.5834

$followUpTime
[1] 12

$analysisTime
         [,1]
[1,] 16.79806
[2,] 29.38334

Staggered patient entry

Assume piecewise accrual, i.e. different enrollment rates over time:

  • 0 to <3 months: 15 patients/month
  • 3 to <6 months: 20 patients/month
  • 6 to 16 months: 25 patients/month

Staggered entry for Case 2 setting

designGS |>
    getSampleSizeSurvival(
        median1 = 18,
        median2 = 12,
        dropoutRate1 = 0.05,
        dropoutRate2 = 0.05,
        dropoutTime = 12,
        accrualTime = c(0, 3, 6, 16),
        accrualIntensity = c(15, 20, 25)
    ) |>
    fetch(
        maxNumberOfSubjects,
        followUpTime,
        analysisTime
    )
$maxNumberOfSubjects
[1] 355

$followUpTime
[1] 23.60973

$analysisTime
         [,1]
[1,] 18.83368
[2,] 39.60973

Staggered entry via list notation

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

designGS |>
    getSampleSizeSurvival(
        median1 = 18,
        median2 = 12,
        dropoutRate1 = 0.05,
        dropoutRate2 = 0.05,
        dropoutTime = 12,
        accrualTime = accrualTimetList
    ) |>
    fetch(
        maxNumberOfSubjects,
        followUpTime,
        analysisTime
    )
$maxNumberOfSubjects
[1] 355

$followUpTime
[1] 23.60973

$analysisTime
         [,1]
[1,] 18.83368
[2,] 39.60973

Staggered entry with maximum subjects

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

designGS |>
    getSampleSizeSurvival(
        median1 = 18,
        median2 = 12,
        dropoutRate1 = 0.05,
        dropoutRate2 = 0.05,
        dropoutTime = 12,
        accrualTime = accrualTimetList,
        maxNumberOfSubjects = 355
    ) |>
    fetch(
        accrualTime,
        followUpTime,
        analysisTime
    )
$accrualTime
[1]  3  6 16

$followUpTime
[1] 23.60973

$analysisTime
         [,1]
[1,] 18.83368
[2,] 39.60973

Staggered entry for Case 3 setting

designGS |>
    getSampleSizeSurvival(
        median1 = 18,
        median2 = 12,
        accrualTime = c(0, 3, 6),
        accrualIntensity = c(15, 20, 25),
        followUpTime = 12
    ) |>
    fetch(
        maxNumberOfSubjects,
        accrualTime,
        analysisTime
    )
$maxNumberOfSubjects
[1] 420.3458

$accrualTime
[1]  3.00000  6.00000 18.61383

$analysisTime
         [,1]
[1,] 18.28469
[2,] 30.61383

Verify results by simulation

Planning assumptions

Assume the trial is planned with:

  • 257 events
  • 400 patients
  • accrual until month 16 at 25/month (with the staggered structure above)

We compare asymptotic and simulation-based results.

Analytical power for fixed event plan

designGS |>
    getPowerSurvival(
        median1 = 18,
        median2 = 12,
        dropoutRate1 = 0.05,
        dropoutRate2 = 0.05,
        dropoutTime = 12,
        accrualTime = c(0, 3, 6, 16),
        accrualIntensity = c(15, 20, 25),
        maxNumberOfEvents = 257,
        directionUpper = FALSE # Important!
    ) |>
    fetch(analysisTime, overallReject)
$analysisTime
         [,1]
[1,] 18.85699
[2,] 39.74904

$overallReject
[1] 0.9005251

Simulation for the same setting

designGS |>
    getSimulationSurvival(
        median1 = 18,
        median2 = 12,
        dropoutRate1 = 0.05,
        dropoutRate2 = 0.05,
        dropoutTime = 12,
        maxNumberOfIterations = 10000,
        accrualTime = c(0, 3, 6, 16),
        accrualIntensity = c(15, 20, 25),
        plannedEvents = c(129, 257), # Note this different input here!
        directionUpper = FALSE,
        seed = 12345
    ) |>
    fetch(analysisTime, overallReject)
$analysisTime
         [,1]
[1,] 18.91939
[2,] 39.63640

$overallReject
[1] 0.9054

The simulation result should closely match the analytic power.

Adaptive survival design with sample size reassessment

Modified scenario

Now assume:

  • Events increase up to 10-fold is allowed
    • Need to make sure that maximum number of subjects is large enough to allow for this
  • Conditional power target is 90%
  • Reassessment is based on observed hazard ratios

For this SSR setting, we use an inverse normal design.

Starting point: Inverse normal design object

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

Design parameters and output of inverse normal combination test design

User defined parameters

  • Type of design: O’Brien & Fleming type alpha spending
  • Maximum number of stages: 2
  • Stages: 1, 2
  • Type II error rate: 0.1000

Derived from user defined parameters

  • Information rates: 0.500, 1.000
  • Futility bounds (non-binding): -Inf

Default parameters

  • Significance level: 0.0250
  • Binding futility: FALSE
  • Test: one-sided
  • Tolerance: 1e-08
  • Type of beta spending: none

Output

  • Cumulative alpha spending: 0.001525, 0.025000
  • Critical values: 2.963, 1.969
  • Stage levels (one-sided): 0.001525, 0.024500

Without SSR: reference result

designIN |>
    getPowerSurvival(
        median1 = 18,
        median2 = 12,
        accrualTime = c(0, 16),
        maxNumberOfSubjects = 2000,
        maxNumberOfEvents = 257,
        directionUpper = FALSE,
        dropoutRate1 = 0.05,
        dropoutRate2 = 0.05,
        dropoutTime = 12
    ) |>
    fetch(analysisTime, overallReject)
$analysisTime
          [,1]
[1,]  6.945121
[2,] 10.081824

$overallReject
[1] 0.9005251

With SSR: simulation result

designIN |>
    getSimulationSurvival(
        median1 = 18,
        median2 = 12,
        maxNumberOfIterations = 10000,
        accrualTime = c(0, 16),
        maxNumberOfSubjects = 2000,
        plannedEvents = c(129, 257),
        directionUpper = FALSE,
        dropoutRate1 = 0.05,
        dropoutRate2 = 0.05,
        dropoutTime = 12,
        conditionalPower = 0.9,
        # no sample size reduction is allowed:
        minNumberOfEventsPerStage = c(NA, 128),
        # up to 10-fold increase:
        maxNumberOfEventsPerStage = 10 * c(NA, 128),
        seed = 23456
    ) |>
    fetch(analysisTime, overallReject)
$analysisTime
          [,1]
[1,]  6.954588
[2,] 15.501731

$overallReject
[1] 0.9855

Type I error control under data-driven SSR

We now compare analysis methods under \(H_0\) (hazardRatio = 1):

  • Inverse normal method (appropriate for data-driven SSR)
  • Group sequential method (inappropriate for this SSR setting, only shown as a contrast here)

Simulation: inverse normal vs group sequential

designIN <- getDesignInverseNormal(
    kMax = 2,
    typeOfDesign = "asOF"
)
designIN |>
    getSimulationSurvival(
        hazardRatio = 1,
        maxNumberOfIterations = 10000,
        accrualTime = c(0, 16),
        maxNumberOfSubjects = 2000,
        plannedEvents = c(129, 257),
        directionUpper = FALSE,
        dropoutRate1 = 0.05,
        dropoutRate2 = 0.05,
        dropoutTime = 12,
        conditionalPower = 0.9,
        minNumberOfEventsPerStage = c(NA, 128),
        maxNumberOfEventsPerStage = 10 * c(NA, 128),
        seed = 34567
    ) |>
    fetch(overallReject)
overallReject 
       0.0222 
designGS |>
    getSimulationSurvival(
        hazardRatio = 1,
        maxNumberOfIterations = 10000,
        accrualTime = c(0, 16),
        maxNumberOfSubjects = 2000,
        plannedEvents = c(129, 257),
        directionUpper = FALSE,
        dropoutRate1 = 0.05,
        dropoutRate2 = 0.05,
        dropoutTime = 12,
        conditionalPower = 0.9,
        minNumberOfEventsPerStage = c(NA, 128),
        maxNumberOfEventsPerStage = 10 * c(NA, 128),
        seed = 45678
    ) |>
    fetch(overallReject)
overallReject 
       0.0369 

Non-proportional hazards

Simulation under non-proportional hazards

  • Assume a delayed treatment effect with a delay of 6 months, i.e. the active active treatment does not affect the hazard in the first 6 months but still reduces it by 1/3 after the delay
  • We can specify this with piecewise exponential survival times with a change point at 6 months:
designIN |> getSimulationSurvival(
    median1 = c(18, 12),
    median2 = 12,
    maxNumberOfIterations = 10000,
    accrualTime = c(0, 16),
    maxNumberOfSubjects = 2000,
    plannedEvents = c(129, 257),
    directionUpper = FALSE,
    conditionalPower = 0.9,
    dropoutRate1 = 0.05,
    dropoutRate2 = 0.05,
    dropoutTime = 12,
    minNumberOfEventsPerStage = c(NA, 128),
    maxNumberOfEventsPerStage = 10 * c(NA, 128),
    seed = 56789
)

Simulation under non-proportional hazards (cont’d)

  • As can be seen from the output above, power drops substantially
  • If one wants to maintain the 90% power target, one needs to increase the maximum number of events
  • Note that in this example the power does not only depend on the number of events, but also on the speed of recruitment!
eventsMax <- 350 # could be further decreased
designIN |>
    getSimulationSurvival(
        piecewiseSurvivalTime = c(0, 6), # new!
        lambda1 = c(log(2) / 12, log(2) / 18), # delay!
        lambda2 = rep(log(2) / 12, 2),
        maxNumberOfIterations = 10000,
        accrualTime = c(0, 16),
        maxNumberOfSubjects = 2000,
        plannedEvents = c(0.6 * eventsMax, eventsMax), # Increased
        directionUpper = FALSE,
        dropoutRate1 = 0.05,
        dropoutRate2 = 0.05,
        dropoutTime = 12,
        conditionalPower = 0.9,
        minNumberOfEventsPerStage = c(NA, 128),
        maxNumberOfEventsPerStage = 10 * c(NA, 128),
        seed = 56789
    ) |>
    fetch(overallReject)
overallReject 
       0.9211 

Futility analyses

Futility analyses

  • Similar as in the continuous endpoint example, we can also consider different futility analysis scenarios here
  • Let’s say we want to vary both the timing (number of events) as well as the futility boundary, specifying it on the conditional power scale for easier interpretation
  • The starting point is the inverse normal design object, where we now in addition specify explicitly the information rate and the futility bound for the interim analysis:
getDesignInverseNormal(
    kMax = 2,
    typeOfDesign = "asOF",
    beta = 0.1,
    informationRates = c(0.5, 1),
    futilityBounds = 0.5,
    futilityBoundsScale = "condPowerAtObserved"
)

Design parameters and output of inverse normal combination test design

User defined parameters

  • Type of design: O’Brien & Fleming type alpha spending
  • Type II error rate: 0.1000
  • Futility bounds (non-binding): 1.392

Derived from user defined parameters

  • Maximum number of stages: 2
  • Stages: 1, 2

Default parameters

  • Information rates: 0.500, 1.000
  • Significance level: 0.0250
  • Binding futility: FALSE
  • Test: one-sided
  • Tolerance: 1e-08
  • Type of beta spending: none

Output

  • Cumulative alpha spending: 0.001525, 0.025000
  • Critical values: 2.963, 1.969
  • Stage levels (one-sided): 0.001525, 0.024500

Futility analyses (cont’d)

Let’s now define the different scenarios:

scenarios <- list(
    list(futilityBounds = 0.6, infRate = 0.6),
    list(futilityBounds = 0.5, infRate = 0.5),
    list(futilityBounds = 0.4, infRate = 0.4),
    list(futilityBounds = 0.3, infRate = 0.4),
    list(futilityBounds = 0.2, infRate = 0.4)
)

Based on this we can calculate calculate required sample sizes and corresponding study durations as well as other trial characteristics for the different scenarios:

results <- scenarios |>
    lapply(function(scenario) {
        getDesignInverseNormal(
            kMax = 2,
            typeOfDesign = "asOF",
            beta = 0.1,
            informationRates = c(scenario$infRate, 1),
            futilityBounds = scenario$futilityBounds,
            futilityBoundsScale = "condPowerAtObserved"
        ) |> getSampleSizeSurvival(
            median1 = 18,
            median2 = 12,
            dropoutRate1 = 0.05,
            dropoutRate2 = 0.05,
            dropoutTime = 12,
            accrualTime = c(0, 16),
            maxNumberOfSubjects = 1000
        )
    })

Futility analyses (cont’d)

Let’s look at the results, for example the maximal study duration, the futility bound on the hazard ratio scale, and the expected number of subjects under \(H_1\):

sapply(results, function(result) {
    result |>
        fetch(
            "Maximal study duration",
            "Futility bounds (treatment effect scale)",
            "Expected number of subjects under H1"
        )
})
                           [,1]      [,2]      [,3]      [,4]      [,5]     
maxStudyDuration           18.19997  18.28056  18.51473  17.59329  16.76924 
futilityBoundsEffectScale  0.7983172 0.8134792 0.8323314 0.8454797 0.8634906
expectedNumberOfSubjectsH1 886.2009  881.4691  912.2011  915.0433  919.939  

Note that you can find out the result statistics names via looking at the print output of a typical result object, e.g. results[[1]].

Key takeaways

  • Different accrual specifications can describe the same recruitment process and be flexibly adapted depending on the desired inputs
  • Analytical formulas and simulation results align closely in this example (and more generally in large enough trials where the asymptotics work well)
  • For data-driven sample size re-estimation, the inverse normal framework maintains Type I error
  • Non-proportional hazards can be accounted for using simulations
  • Futility analyses can be flexibly defined and evaluated in terms of key trial characteristics