Outcome and dropout distributions

We begin with information about the median time-to-event in thecontrol group, dropout rate, hazard ratios under the null and alternatehypotheses for experimental therapy compared to control, and the desiredType I and II error rates.

# Median control time-to-eventmedian<-12# Exponential dropout rate per unit of timeeta<- .001# Hypothesized experimental/control hazard ratio# (alternate hypothesis)hr<- .75# Null hazard ratio (1 for superiority, >1 for non-inferiority)hr0<-1# Type I error (1-sided)alpha<- .025# Type II error (1-power)beta<- .1

Enrollment and trial duration

Next, we plan the trial duration and the enrollment pattern. Thereare two basic methods for doing this. TheLachinand Foulkes (1986) method demonstrated here fixes the enrollmentpattern and duration as well as the trial duration and changes theabsolute enrollment rates to obtain desired power. The alternaterecommended method is along the lines ofKim andTsiatis (1990), fixing the enrollment rates and follow-upduration, varying the total trial duration to power the design; thiswill also be demonstrated below.

# Study durationT<-36# Follow-up duration of last patient enrolledminfup<-12# Enrollment period durationsR<-c(1,2,3,4)# Relative enrollment rates during above periodsgamma<-c(1,1.5,2.5,4)# Randomization ratio, experimental/controlratio<-1

Deriving design with no interim analyses

The above information is sufficient to design a trial with no interimanalyses. Note that when callingnSurv(), we transform themedian time-to-event (\(m\)) to anexponential event rate (\(\lambda\))with the formula\[\lambda=\log(2)/m.\]

library(gsDesign)x<-nSurv(R = R,gamma = gamma,eta = eta,minfup = minfup,T = T,lambdaC =log(2)/ median,hr = hr,hr0 = hr0,beta = beta,alpha = alpha)

A textual summary of this design is given by printing it. For thegroup sequential design shown later, much more complete formatted outputwill be shown.

x#> Fixed design, two-arm trial with time-to-event#> outcome (Lachin and Foulkes, 1986).#> Solving for:  Accrual rate#> Hazard ratio                  H1/H0=0.75/1#> Study duration:                   T=36#> Accrual duration:                   24#> Min. end-of-study follow-up: minfup=12#> Expected events (total, H1):        507.1519#> Expected sample size (total):       775.0306#> Accrual rates:#>      Stratum 1#> 0-1     9.2818#> 1-3    13.9227#> 3-6    23.2045#> 6-24   37.1272#> Control event rates (H1):#>       Stratum 1#> 0-Inf    0.0578#> Censoring rates:#>       Stratum 1#> 0-Inf     0.001#> Power:                 100*(1-beta)=90%#> Type I error (1-sided):   100*alpha=2.5%#> Equal randomization:          ratio=1

Varying enrollment duration to power trial

If we had setT = NULL above, the specified enrollmentrates would not be changed but enrollment duration would be adjusted toachieve desired power. For the low enrollment rates specified ingamma above, this would have resulted in a long trial.

# THIS CODE IS EXAMPLE ONLY; NOT EXECUTED HEREnSurv(R = R,gamma = gamma,eta = eta,minfup = minfup,T =NULL,# This was changedlambdaC =log(2)/ median,hr = hr,hr0 = hr0,beta = beta,alpha = alpha)

Additional parameters

All the parameters above are used. We set up the number of analyses,timing and spending function parameters. These deserve careful attentionfor every trial and tend to be somewhat customized to be fit-for-purposeaccording to all those involved in designing the trial. Here the choicesconsidered the following:

• There is an early look intended primarily for futility
  • There is no desire to stop this early for lack of a positivetrend
  • An earlier interim might be too soon to stop for a meaningful safetyproblem
  • There is no desire to stop for efficacy, so the extreme effect sizeand p-value required to cross the bound is considered appropriate
• The second interim is primarily for efficacy
  • The data will be relatively mature, but a meaningful length of timebefore planned final analysis
  • The treatment effect at the efficacy bound will likely be at leastthe treatment effect used to power the trial
  • The Lan-DeMets spending function approximating an O’Brien-FlemingboundLan and DeMets (1983) is oftenacceptable to regulators for early stopping

# Number of analyses (interim + final)k<-3# Timing of interim analyses (k-1 increasing numbers >0 and <1).# Proportion of final events at each interim.timing<-c(.25, .75)# Efficacy bound spending function.# We use Lan-DeMets spending function approximating O'Brien-Fleming bound.# No parameter required for this spending function.sfu<- sfLDOFsfupar<-NULL# Futility bound spending functionsfl<- sfHSD# Futility bound spending parameter specificationsflpar<--7

Type II error (1-power) may be set up differently than for a fixeddesign so that more meaningful futility analyses can be performed duringthe course of the trial.

# Type II error = 1 - Powerbeta<- .15

Generating the design

Now we are prepared to generate the design.

# Generate designx<-gsSurv(k = k,timing = timing,R = R,gamma = gamma,eta = eta,minfup = minfup,T = T,lambdaC =log(2)/ median,hr = hr,hr0 = hr0,beta = beta,alpha = alpha,sfu = sfu,sfupar = sfupar,sfl = sfl,sflpar = sflpar)

Textual summary

The design summary is:

cat(summary(x))

Asymmetric two-sided group sequential design with non-bindingfutility bound, 3 analyses, time-to-event outcome with sample size 676and 443 events required, 85 percent power, 2.5 percent (1-sided) Type Ierror to detect a hazard ratio of 0.75. Enrollment and total studydurations are assumed to be 24 and 36 months, respectively. Efficacybounds derived using a Lan-DeMets O’Brien-Fleming approximation spendingfunction (no parameters). Futility bounds derived using aHwang-Shih-DeCani spending function with gamma = -7.

An important addition not provided above is that the mediantime-to-event is assumed to be 12 months in the control group.

Tabular summaries

Following are the enrollment rates required to power the trial.

library(gt)library(tibble)tibble(Period =paste("Month",rownames(x$gamma)),Rate =as.numeric(x$gamma))%>%gt()%>%tab_header(title ="Enrollment rate requirements")

Next we provide a tabular summary of bounds for the design. We haveadded extensive footnoting to the table, which may or may not berequired for your design. However, as seen here it makes many choicesfor design parameters and properties transparent. No attempt has beenmade to automate this, but it may be worth considering for a template ifyou wish to make the same choice across many trials. Note that theexclude argument forgsBoundSummary() allowsadditional descriptions for design bounds such as conditional orpredictive power; see the help file for details or just provideexclude = NULL togsBoundSummary() to see alloptions.

# Footnote text for tablefootnote1<-"P{Cross} is the probability of crossing the given bound (efficacy or futility) at or before the given analysis under the assumed hazard ratio (HR)."footnote2<-" Design assumes futility bound is discretionary (non-binding); upper boundary crossing probabilities shown here assume trial stops at first boundary crossed and thus total less than the design Type I error."footnoteHR<-"HR presented is not a requirement, but an estimate of approximately what HR would be required to cross each bound."footnoteM<-"Month is approximated given enrollment and event rate assumptions under alternate hypothesis."# Spending function footnotesfootnoteUS<-"Efficacy bound set using Lan-DeMets spending function approximating an O'Brien-Fleming bound."footnoteLS<-paste("Futility bound set using ", x$lower$name," beta-spending function with ",  x$lower$parname,"=", x$lower$param,".",sep ="")# Caption text for tablecaption<-paste("Overall survival trial design with HR=", hr,", ",100* (1- beta),"% power and ",100* alpha,"% Type I error",sep ="")

gsBoundSummary(x)%>%gt()%>%tab_header(title ="Time-to-event group sequential design")%>%cols_align("left")%>%tab_footnote(footnoteUS,locations =cells_column_labels(columns =3))%>%tab_footnote(footnoteLS,locations =cells_column_labels(columns =4))%>%tab_footnote(footnoteHR,locations =cells_body(columns =2,rows =c(3,8,13)))%>%tab_footnote(footnoteM,locations =cells_body(columns =1,rows =c(4,9,14)))%>%tab_footnote(footnote1,locations =cells_body(columns =2,rows =c(4,5,9,10,14,15)))%>%tab_footnote(footnote2,locations =cells_body(columns =2,rows =c(4,9,14)))

Summary plots

Several plots are available to summarize a design; see help forplot.gsDesign(); one easy way to see how to generate eachis by checking plots and code generated by theShiny interface. Thepower plot is information-rich, but also requires some explanation;thus, we demonstrate here.

The solid black line represents the trial power by effect size. Powerat interim 1 is represented by the black dotted line. Cumulative powerat interim 2 is represented by the black dashed line. The red dottedline is 1 minus the probability of crossing the futility bound on thepercentage scale. The red dashed line is 1 minus the cumulativeprobability of crossing the futility bound by interim 2.

library(ggplot2)library(scales)plot(x,plottype ="power",cex = .8,xlab ="HR")+scale_y_continuous(labels = scales::percent)

# Number of events (final is still planned number)n.I<-c(115,364,ceiling(x$n.I[x$k]))

xu<-gsDesign(alpha = x$alpha,beta = x$beta,test.type = x$test.type,maxn.IPlan = x$n.I[x$k],n.I = n.I,sfu = sfu,sfupar = sfupar,sfl = sfl,sflpar = sflpar,delta = x$delta,delta1 = x$delta1,delta0 = x$delta0)

gsBoundSummary(  xu,deltaname ="HR",logdelta =TRUE,Nname ="Events",exclude =c("Spending","B-value","CP","CP H1","PP","P(Cross) if HR=1","P(Cross) if HR=0.75"  ))%>%gt()%>%cols_align("left")%>%tab_header(title ="Time-to-event group sequential bound guidance",subtitle ="Bounds updated based on event counts through IA2"  )%>%tab_footnote("Nominal p-value required to establish statistical significance.",locations =cells_body(columns =3,rows =c(2,5,8))  )%>%tab_footnote("Interim futility guidance based on observed HR is non-binding.",locations =cells_body(columns =4,rows =c(3,6))  )%>%tab_footnote("HR bounds are approximations; decisions on crossing are based solely on p-values.",locations =cells_body(column =2,rows =c(3,6,9))  )

Evaluating interim results

We recommend 3 things to present to summarize results in addition tostandard summaries such as the logrank p-value, hazard ratio based onthe Cox model, median time-to-event and Kaplan-Meier curves by treatmentgroup.

• Conditional power at the current trend, under the null hypothesis ofno treatment difference (conditional error), and under the alternativehypothesis.
• Predictive power, which is conditional power averaged over aposterior distribution for the treatment effect.
• A B-value plot to evaluate the trend in test statistics towards apositive conclusion.

For these summaries, we will assume the updated interim event countsused above along with interim Z-values of 0.25 and 2 at interim 1 andinterim 2, respectively.

Z<-c(0.25,2)

Conditional power at interim analysis 2 is computed for the currenttrend, under the null hypothesis (HR=1), and under the alternatehypothesis (HR=0.75 in this case) as follows:

gsCP(x = xu,# Updated designi =2,# Interim analysis 2zi = Z[2]# Observed Z-value for testing)$upper$prob#>           [,1]     [,2]      [,3]#> [1,] 0.6599398 0.301728 0.7764629

Predictive power incorporates uncertainty into the above conditionalpower evaluation. The computation assumes a prior distribution for thetreatment effect and then updates to a posterior distribution for thetreatment effect based on the most recent interim result. Theconditional probability of a positive finding is then averaged accordingto this posterior. We specify a normal prior for the standardized effectsize using thegsDesign::normalGrid() function. We select aweak prior with mean half-way between the alternative(x$delta) and null (0) hypotheses and a variance equivalentto observing 5% (=1/20) of the targeted events at the final analysis;the following shows that the standard deviation for the prior is wellover twice the mean, so the prior is relatively weak.

prior<-normalGrid(mu = x$delta/2,sigma =sqrt(20/max(x$n.I)))cat(paste(" Prior mean:",round(x$delta/2,3),"\n Prior standard deviation",round(sqrt(20/ x$n.fix),3),"\n"))#>  Prior mean: 0.072#>  Prior standard deviation 0.215

Now based on the interim 2 result, we compute the predictive power ofa positive final analysis.

gsPP(x = xu,# Updated designi =2,# Interim analysis 2zi = Z[2],# Observed Z-value for testingtheta = prior$z,# Grid points for above priorwgts = prior$wgts# Weights for averaging over grid)#> [1] 0.6407376

A B-value (Proschan, Lan, and Wittes(2006)) is a Z-value multiplied by the square root of theinformation fraction (interim information divided by final plannedinformation. In the plot below on the B-value scale, we present theefficacy bounds at each analysis in black, futility guidance in red, theobserved interim tests in blue connected by solid lines, and a dashedblue line to project the final result. Under a constant treatment effect(proportional hazards for a time-to-event outcome tested with a logranktest) the blue line behaves like observations from a Brownian motionwith a linear trend (“constant drift”). While a comparable Z-value plotwould have the effect increasing with the square root of the number ofevents, the B-value plot trend is linear in the event count. The trendis proportional to the logarithm of the underlying hazard ratio. Theprojected final test is based on the dashed line which represents alinear trend based on the most recent B-value computed; this projectionis what was used in the conditional power calculation under the currenttrend that was computed above.

maxx<-450# Max for x-axis specified by userylim<-c(-1,3)# User-specified y-axis limitsanalysis<-2# Current analysis specified by user# Following code should require no further changesplot(  xu,plottype ="B",base =TRUE,xlim =c(0, maxx),ylim = ylim,main ="B-value projection",lty =1,col =1:2,xlab ="Events")N<-c(0, xu$n.I[1:analysis])B<-c(0, Z*sqrt(xu$timing[1:analysis]))points(x = N,y = B,col =4)lines(x = N,y = B,col =4)slope<- B[analysis+1]/ N[analysis+1]Nvals<-c(N[analysis+1],max(xu$n.I))lines(x = Nvals,y = B[analysis+1]+c(0, slope* (Nvals[2]- Nvals[1])),col =4,lty =2)