getSimulationSurvival, thetaH0 and conditional power

[grlm1]

Hi! First of all thank you for this very nice and needed package!

I am trying to simulate an adaptive design using the function getSimulationSurvival(). I would also like to use thetaH0=0.9 instead of the default thetaH0=1. I am wondering if this parameter is respected during the simulation and/or in the results.

library(rpact)
library(dplyr)
dIN <- getDesignInverseNormal(kMax=2,
                              alpha=0.025,
                              beta=0.1,
                              sided=1,
                              typeOfDesign = "noEarlyEfficacy",
                              informationRates = c(19/125, 1),
                              futilityBounds = c(-1))

dIN_sim <-  getSimulationSurvival(design = dIN, 
                                  pi2 = 0.003,
                                  eventTime = 12,
                                  hazardRatio = c(0.5),
                                  dropoutRate1 = 0.15, dropoutRate2 = 0.15, dropoutTime = 12,
                                  accrualTime = c(0, 20),
                                  accrualIntensity = c(1000),
                                  thetaH0 = 0.9,
                                  directionUpper = F,
                                  plannedEvents = c(19, 125),              #cumulated
                                  maxNumberOfEventsPerStage = c(NA, 200),  #max events
                                  minNumberOfEventsPerStage = c(NA, 106),  #min events (per stage, not cumulative!)
                                  conditionalPower = 0.9,
                                  seed = 1910)

Since I would like to adapt the simulation manually (specifically, restart subject accrual at the timing of the interim based on conditional power and/or pi2), I tried to replicate scenarios from the simulation data using the smaller functions getStageResults() and getConditionalPower() .

simulation_replicates <- getData(dIN_sim) %>% mutate(iterationNumber = factor(iterationNumber))
interesting_case <- simulation_replicates %>% filter(iterationNumber == 12)
interesting_case  #increased sample size slightly from 106 to 131, i.e. overall 15o cumulated events in stage 2

which gives me

iterationNumber stageNumber         pi1   pi2 hazardRatio analysisTime numberOfSubjects overallEvents1 overallEvents2 eventsPerStage rejectPerStage
1              12           1 0.001501127 0.003         0.5     15.34246            15342              7             12             19              0
2              12           2 0.001501127 0.003         0.5     84.62074            20000             49            101            131              1
  eventsNotAchieved futilityPerStage testStatistic logRankStatistic conditionalPowerAchieved pValuesSeparate trialStop hazardRatioEstimateLR
1                 0                0     0.9191049        0.9191049                       NA    0.1790203579     FALSE             0.6559214
2                 0                0     3.6488241        3.6648582                0.7509866    0.0001774007      TRUE             0.5496526

Plugging in the results of the getSimulationSurvival() interim analysis result into getDataset() followed by getStageResults() and getConditionalPower() to the best of my ability:

#with thetah0=0.9
cp_with_thetah0_09 <- getConditionalPower(stageResults = getStageResults(dIN, dataInput = getDataset(cumulativeEvents=c(19), 
                                                                               cumulativeLogRanks = c(-0.9191049),  # minus! 
                                                                               cumulativeAllocationRatios = c(1)), 
                                                   directionUpper = F, 
                                                   stage = 1,
                                                   thetaH0 = 0.9),
                    thetaH1 = 0.6559214,   #use HR of stage 1
                    nPlanned = 131) #adaptively determined ssr for stage 2

#with thetah0=1
cp_with_thetah0_1 <- getConditionalPower(stageResults = getStageResults(dIN, dataInput = getDataset(cumulativeEvents=c(19), 
                                                                              cumulativeLogRanks = c(-0.9191049),  # minus! 
                                                                              cumulativeAllocationRatios = c(1)), 
                                                                         directionUpper = F, 
                                                                         stage = 1,
                                                                         thetaH0 = 1),
                                          thetaH1 = 0.6559214,   #use HR of stage 1
                                          nPlanned = 131)      #adaptively determined ssr for stage 2

Comparing the results for reported conditional power:

> cp_with_thetah0_09$conditionalPower
[1]        NA 0.4896989
> cp_with_thetah0_1$conditionalPower
[1]        NA 0.7497147
> interesting_case$conditionalPowerAchieved
[1]        NA 0.7509866

which made me suspect that the reported conditional power of getSimulationSurvival() does not respect the thetaH0 parameter. Or is $conditionalPower not to be confused with $conditionalPowerAchieved? Additionally, why do I need to plug in the negative LR statistic into the getStageResults(getDataset()) function? I suspect the reason is the directionUpper argument, but I fail to see why the getSimulationSurvival() does not report the negative LR's that I would expect for HR's < 1.

Thank you!

[fpahlke]

Hi Tim,

Thank you for your interest in rpact and pointing out this issue. In fact, you find a bug and the non-inferiority bound is not taken into account in the sample size reassessment procedure (at least, the output is not correct). The bug will be fixed as quickly as possible.

Concerning your second question: the test statistic in the simulation output is adjusted according to the direction of the alternative. If you want to use it in getStageResults you should either withdraw the adjustment or set directionUpper = TRUE.

Hope this helps, we will inform you when the correction will be available on GitHub.

Kind regards, Gernot and Friedrich

[fpahlke]

The solution is now available in the development branch dev/3.5.0. The next release on CRAN is planned for the end of January 2024.