---
  title: "RMarkdown version of lab CAST"
  author: "Constantin T Yiannoutsos"
  date: "March 21, 2018"
  output: html_document
---
  
```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
```

```{r, packages, echo=FALSE, results="hide", message=FALSE}
options(width=72)
```

# Introduction: The R package ```gsDesign```
First, we import the library ```gsDesign``` into R.
```{r, enrollment, message=FALSE}
################################################
### Sequential monitoring of clinical trials ###
################################################
library(gsDesign)
```
# Defining user-specified spending functions 

The function ```gsDesign``` can carry out lots of analyses. Here we focus on the definition of a user-specified spending fuction.

The CAST study had a function that spent half the alpha level (of 2.5%) during the interim analyses prior to the final analysis, and half of the alpha level during the final analysis.

We first define the spending function:
```{r, sf}
sfcast <- function(alpha,  t,  param)
{  

  # set spending using constant alpha/2 until last interim then alpha/2 at final
  x <- list(name="CAST example", param=param, parname=NULL, 
            sf=sfcast, spend=cumsum(c(rep(alpha/(2*(length(t)-1)),length(t)-1), alpha/2)))  
  
  class(x) <- "spendfn"
  
  x
}
```
Notice that the object that is being created ($x$) has a number of attributes, such as ```name```, ```param``` and it must be of class ```spendfn``` to be recognized by the package.

Now let's try it out!

First we will plot the spending fuunction to see hwo the alpha level is being spent

```{r, sfexample,  fig.width=4.5, fig.height=4.5, message=FALSE, fig.align='center', fig.cap='<b>Figure 1. </b>Alpha spending function based on user-specified spending function, assuming  &alpha;=0.025.'}
t<-0:20/20
plot(t, sfcast(0.025, t, NULL)$spend, type="l", 
     xlab="Proportion of information", 
     ylab="Cumulative proportion of total spending")
```
# Using ```sfcast``` in the CAST trial example

Now, let's define the CAST trial bounds:
```{r, CASTbounds, fig.width=4.5, fig.height=4.5, message=FALSE, fig.align='center', fig.cap='<b>Figure 2. </b> Group sequential bounds in the CAST study.'}
x <- gsDesign(k=length(t), test.type=2, sfu=sfcast, alpha=0.025)
plot(x$timing, x$upper$bound, type="l", ylim=c(-3.5,3.5) ,
     ylab="Z(t)", xlab="Trial fraction")
lines(x$timing, x$lower$bound)
```
Now recall that the first DSMB review of the CAST study was after 32/425 events were observed (i.e., $t=0.075$) and the log-rank test statistic was $Z(t)=-2.82$. During the second DSMB meeting, 42/300 events were observed ($t=14\%$) after the expected number of total events was revised to 300 from 425.  During that time, $Z(t)=-3.20$.

To see whether we have crossed a boundary in either case, we simply plot these two points on the previous plot.

```{r, CASTreviews, fig.width=4.5, fig.height=4.5, message=FALSE, fig.align='center', fig.cap='<b>Figure 3. </b> Group sequential bounds in the CAST study with log-rank statistics at the two DSMB reviews.'}
x <- gsDesign(k=length(t), test.type=2, sfu=sfcast, alpha=0.025)
plot(x$timing, x$upper$bound, type="l", ylim=c(-3.5,3.5) ,
     ylab="Z(t)", xlab="Trial fraction")
lines(x$timing, x$lower$bound)
text(32/425, -2.82, "X")
text(42/300, -3.2, "X", col="red")
```