##
#####
#####	Solutions to Session 3 Problems
#####
##### Author: D. Gillen
#####
##
#####
#####	Problem 1: Analytic results vs. simulation
#####
##
library( RCTdesign )
##
#####	Pocock Design
##
poc <- seqDesign( 	prob.model="normal", arms=1, null=0, alt=1, var=2, 
					nbr.analyses=4, P=0.5, test.type="two.sided", 
					early.stopping="alternative", sample.size=100, power="calculate" )
poc
seqBoundary( poc, scale="Z" )	## Simulated critical value was 2.635 (z-scale)


##
#####	OBF Design
##
obf <- update( poc, P=1 )
obf
seqBoundary( poc, scale="S" )	## Simulated critical value was 28.643 (partial sum scale)

##
#####	Plot of stopping boundaries on Z-scale
##
plot( poc, obf, display.scale="Z" )

##
#####	Comparison of operating characteristics (ASN and power)
seqOC( poc, theta=(0:10)/10 )
seqOC( poc, theta=(1:10)/10 )$power.lower + seqOC( obf, theta=(1:10)/10 )$power.upper

seqOC( obf, theta=(0:10)/10 )
seqOC( obf, theta=(1:10)/10 )$power.lower + seqOC( obf, theta=(1:10)/10 )$power.upper

##
#####	Results from simulation
##
#> cbind( mu, pocPower, obfPower, fixedPower )
#       mu pocPower obfPower fixedPower
# [1,] 0.1  0.09170  0.10673    0.10954
# [2,] 0.2  0.23423  0.28665    0.29606
# [3,] 0.3  0.47155  0.55160    0.56329
# [4,] 0.4  0.72856  0.79828    0.80733
# [5,] 0.5  0.90240  0.93933    0.94356
# [6,] 0.6  0.97590  0.98727    0.98862
# [7,] 0.7  0.99622  0.99827    0.99852
# [8,] 0.8  0.99957  0.99983    0.99985
# [9,] 0.9  0.99997  0.99999    0.99999
#[10,] 1.0  1.00000  1.00000    1.00000

#> cbind( mu, pocASN, obfASN, fixedASN )
#       mu pocASN obfASN fixedASN
# [1,] 0.1 96.370 98.621      100
# [2,] 0.2 91.621 95.730      100
# [3,] 0.3 83.165 89.996      100
# [4,] 0.4 71.288 81.206      100
# [5,] 0.5 58.864 71.600      100
# [6,] 0.6 47.965 62.725      100
# [7,] 0.7 39.928 55.765      100
# [8,] 0.8 34.330 50.575      100
# [9,] 0.9 30.552 46.473      100
#[10,] 1.0 28.088 42.887      100


##
#####
#####	Problem 2: Exploring group sequential survival designs
#####
##
#####	Part a
##
#####	Specification of designs
Fixed <- seqDesign( prob.model = "hazard", arms = 2, null.hypothesis = 1., 
							alt.hypothesis = 0.67, ratio = c(1., 1.), nbr.analyses = 1, 
							test.type = "less", power = 0.80, alpha = 0.025 )
Eff11.Fut8 <- update( Fixed, nbr.analyses = 4, P=c(1.1,.8), sample.size=196, power="calculate" )
Eff11.Fut8.4 <- Eff11.Fut8
Eff11.Fut8.2 <- update( Eff11.Fut8.4, nbr.analyses=2 )
Eff11.Fut8.3 <- update( Eff11.Fut8.4, nbr.analyses=3 )
Eff11.Fut8.5 <- update( Eff11.Fut8.4, nbr.analyses=5 )

##
#####	Comparison of power
##
seqPlotPower( Eff11.Fut8.2, Eff11.Fut8.3, Eff11.Fut8.4, Eff11.Fut8.5 )
seqPlotPower( Eff11.Fut8.2, Eff11.Fut8.3, Eff11.Fut8.4, Eff11.Fut8.5, reference=TRUE )

##
#####	Comparison of sample size distribution
##
seqPlotASN( Eff11.Fut8.2, Eff11.Fut8.3, Eff11.Fut8.4, Eff11.Fut8.5 )


##
#####	Part b
##
seqBoundary( Eff11.Fut8.2, scale="E" )
seqBoundary( Eff11.Fut8.3, scale="E" )
seqBoundary( Eff11.Fut8.4, scale="E" )
seqBoundary( Eff11.Fut8.5, scale="E" )

##
#####	Part c
##
Eff11.Fut8.early <- update( Eff11.Fut8.4, sample.size=196*c(.1,.2,.3,1) )
Eff11.Fut8.late <- update( Eff11.Fut8.4, sample.size=196*c(.7,.8,.9,1) )
Eff11.Fut8.early
Eff11.Fut8.late

##
#####	Comparison of power
##
seqPlotPower( Eff11.Fut8.4, Eff11.Fut8.early, Eff11.Fut8.late )
seqPlotPower( Eff11.Fut8.4, Eff11.Fut8.early, Eff11.Fut8.late, reference=TRUE )

##
#####	Comparison of sample size distribution
##
seqPlotASN( Eff11.Fut8.4, Eff11.Fut8.early, Eff11.Fut8.late )

##
##### The above orderings depend upon the boundary shapes!
##
Poc.equal <- update( Fixed, nbr.analyses = 4, P=c(.5,.5), sample.size=196, power="calculate" )
Poc.early <- update( Poc.equal, sample.size=196*c(.1,.2,.3,1) )
Poc.late <- update( Poc.equal, sample.size=196*c(.7,.8,.9,1) )
seqPlotPower(Poc.equal, Poc.early, Poc.late, reference=TRUE)
seqPlotASN( Poc.equal, Poc.early, Poc.late )

##
#####	Part d
##
Eff11.Fut8.Ratio <- update( Eff11.Fut8.4, ratio=2, power=0.7804 )
Eff11.Fut8
Eff11.Fut8.Ratio
seqPHSubjects( Eff11.Fut8.Ratio, controlMedian = 0.75, accrualTime = 3, followupTime = 1 )
ceiling( 278.37 *c(2/3,1/3) )

##
#####
#####	Problem 3: Evaluating designs for the sepsis trial
#####
##
#####
#####	See Emerson, Kittelson, and Gillen, Statistics in Medicine, 2007 
#####	(http://rctdesign.org/TechReports/SIMfreqeval.pdf) for a full 
#####	tutorial evaluating designs for the sepsis trial
#####
##